ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/px/cut.c
Revision: 1.1
Committed: Thu Feb 2 10:49:12 1989 UTC (35 years, 3 months ago) by greg
Content type: text/plain
Branch: MAIN
Log Message:
Initial revision

File Contents

# User Rev Content
1 greg 1.1 /*
2    
3     #ifndef lint
4     static char SCCSid[] = "$SunId$ LBL";
5     #endif
6     CUT - Median bisection (k-d tree) algorithm for color image quantization
7    
8     flaw: doesn't always create as many colors as requested
9     because the binary subdivision is always done in a balanced tree and when a box
10     contains only one color it's unsplittable
11    
12     Paul Heckbert
13     */
14    
15     #include "ciq.h"
16    
17     struct index {short *o; int freq;};
18     /* offset and pixel count for a box */
19    
20     struct pail {struct plum *first; int freak;};
21     /* each pail is a bucket containing plums */
22    
23     struct plum {struct plum *next; int code;};
24     /* a color in a pail's linked list */
25    
26     static short off[len];
27     /* color codes: offsets of colors in hist array */
28    
29    
30     makecm(nw,na) /* subdivide colorspace to generate a colormap*/
31     int nw,*na; /* number of colors wanted, actual */
32     {
33     int i,m,n,freq,weight,sr,sg,sb;
34     short *o;
35     struct index index[257]; /* pointers into off */
36     /* index[i].o is pointer to code for first color in box i */
37    
38     for (o=off, *na=i=0; i<len; i++)
39     if (hist[i]) {
40     (*na)++;
41     *o++ = i; /* initialize offsets */
42     }
43     index[0].o = off; /* start with one big box which contains all *na colors */
44     index[0].freq = xmax*ymax; /* and (sum of hist[i]) pixels */
45     index[1].o = off+*na;
46    
47     /* breadth-first binary subdivision: try to make nw boxes */
48     for (m=1; m<nw; m<<=1) {
49     /*printf("%d ",m);*/
50     index[m<<1].o = index[m].o;
51     for (i=m-1; i>=0; --i) splitbox(index+i,index+(i<<1));
52     }
53     for (n=i=0; i<m && n<nw; i++) if (index[i+1].o>index[i].o) {
54     sr = sg = sb = (freq = index[i].freq)>>1;
55     for (o=index[i].o; o<index[i+1].o; o++) {
56     weight = hist[*o];
57     /* find weighted average of colors in box i */
58     sr += red(*o)*weight;
59     sg += gre(*o)*weight;
60     sb += blu(*o)*weight;
61     }
62     color[0][n] = sr/freq;
63     color[1][n] = sg/freq;
64     color[2][n] = sb/freq;
65     /* printf("color[%d]=(%3d,%3d,%3d) freq=%d\n",
66     n,color[0][n],color[1][n],color[2][n],freq); */
67     n++;
68     }
69     /*printf("[%d empty boxes] from %d to %d colors\n",nw-n,*na,n);*/
70     return n;
71     }
72    
73     splitbox(ii,io) /* split a box in two: half of the pixels from */
74     struct index *ii,*io; /* box ii will go into io, the other half into io+1 */
75     {
76     register short *o,*o1,*o2;
77     register int shift,count,k,freq,r,g,b,r1,g1,b1,r2,g2,b2;
78     struct pail pail[32],*p,*q; /* buckets for sorting colors */
79     struct plum plum[len],*pl=plum,*u; /* colors */
80    
81     o1 = ii[0].o;
82     o2 = ii[1].o;
83     freq = ii[0].freq; /* number of pixels with color in this box */
84     r1 = g1 = b1 = 256; r2 = g2 = b2 = -1;
85     for (o=o1; o<o2; o++) {
86     r = red(*o);
87     g = gre(*o);
88     b = blu(*o);
89     if (r<r1) r1 = r; if (r>r2) r2 = r; /* find min&max r, g, and b */
90     if (g<g1) g1 = g; if (g>g2) g2 = g;
91     if (b<b1) b1 = b; if (b>b2) b2 = b;
92     }
93    
94     /* adaptive partitioning: find dominant dimension */
95     shift = r2-r1>g2-g1?r2-r1>=b2-b1?10:0:g2-g1>=b2-b1?5:0;
96     /* printf("splitbox %3d-%3d %3dx%3dx%3d %c ",
97     o1-off,o2-off,r2-r1,g2-g1,b2-b1,shift==10?'R':shift==5?'G':'B'); */
98    
99     for (p=pail, k=0; k<32; k++, p++) { /* start with empty pails */
100     p->first = 0;
101     p->freak = 0;
102     }
103     for (o=o1; o<o2; o++) {
104     p = pail+(*o>>shift&31); /* find appropriate pail */
105     pl->code = *o;
106     pl->next = p->first;
107     p->first = pl++; /* add new plum to pail */
108     p->freak += hist[*o];
109     }
110    
111     /* find median along dominant dimension */
112     for (count=freq>>1, p=pail; count>0; p++) count -= p->freak;
113    
114     if (p>pail && count+(p-1)->freak<-count) /* back up one? */
115     count += (--p)->freak;
116    
117     io[1].o = o1; /* in case of empty box */
118     for (o=o1, q=pail; o<o2;) { /* dump the pails to reorder colors in box */
119     for (u=q->first; u; u=u->next) *o++ = u->code;
120     if (++q==p) io[1].o = o; /* place partition at offset o */
121     }
122     io[0].o = o1;
123     io[0].freq = (freq>>1)-count;
124     io[1].freq = (freq+1>>1)+count;
125     /* printf(" at %3d %d=%d+%d\n",io[1].o-off,freq,io[0].freq,io[1].freq); */
126     }