ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/cal/cnt.c
(Generate patch)

Comparing ray/src/cal/cnt.c (file contents):
Revision 1.3 by greg, Mon Apr 11 18:08:19 2022 UTC vs.
Revision 1.7 by greg, Fri Apr 22 15:50:34 2022 UTC

# Line 13 | Line 13 | static const char      RCSid[] = "$Id$";
13   #include  <time.h>
14   #include  "random.h"
15  
16 + #ifndef uint16
17 + #define uint16  unsigned short          /* 16-bit unsigned integer */
18 + #endif
19 + #undef uby8
20 + #define uby8  unsigned char             /* 8-bit unsigned integer */
21 +
22   #define MAXDIM  50
23  
24 < char  buf[256];
24 > #define NLEVELS 9                       /* number of tree levels */
25 > #define BRORDER 6                       /* branches/level */
26  
27 + /* Tree branch structure for quick occupancy search */
28 + /* with 9 levels & 6 branches per level, we can store 1.94 Gbits in 259 MBytes (4.5% overhead) */
29 + const struct {
30 +        long    capacity;               /* slots/branch this level */
31 +        long    skip_bytes;             /* bytes until next branch */
32 +        int     cntr_siz;               /* occupancy counter size */
33 + } tree_br[NLEVELS] = {
34 +        {248L, 32L, 1},
35 +        {248L*6, 32L*6+2, 2},
36 +        {248L*6*6, (32L*6+2)*6+2, 2},
37 +        {248L*6*6*6, ((32L*6+2)*6+2)*6+2, 2},
38 +        {248L*6*6*6*6, (((32L*6+2)*6+2)*6+2)*6+3, 3},
39 +        {248L*6*6*6*6*6, ((((32L*6+2)*6+2)*6+2)*6+3)*6+3, 3},
40 +        {248L*6*6*6*6*6*6, (((((32L*6+2)*6+2)*6+2)*6+3)*6+3)*6+3, 3},
41 +        {248L*6*6*6*6*6*6*6, ((((((32L*6+2)*6+2)*6+2)*6+3)*6+3)*6+3)*6+4, 4},
42 +        {248L*6*6*6*6*6*6*6*6, (((((((32L*6+2)*6+2)*6+2)*6+3)*6+3)*6+3)*6+4)*6+4, 4},
43 + };
44  
45 + char  buf[256];                         /* buffer for ordered array output */
46 +
47 +
48   /* Encode integer in string and return pointer to end */
49   static char *
50 < tack(
24 < char  *b,
25 < int  i
26 < )
50 > tack(char *b, long i)
51   {
52          char  *cp;
53          char  *res;
# Line 34 | Line 58 | int  i
58                  *cp++ = '0';
59          else
60                  do {
61 <                        *cp++ = i%10 + '0';
62 <                        i /= 10;
61 >                        *cp++ = i%10L + '0';
62 >                        i /= 10L;
63                  } while (i);
64          res = cp--;
65   #define c i
# Line 51 | Line 75 | int  i
75  
76   /* Loop over dimensions, spitting out buffer after each increment */
77   static void
78 < loop(
55 < int  *n,
56 < char  *b
57 < )
78 > loop(long *n, char *b)
79   {
80 <        int  i;
80 >        long  i;
81  
82          if (n[0] == 0) {
83                  *b = '\0';
# Line 68 | Line 89 | char  *b
89   }
90  
91  
92 + /* Print out shuffled value */
93 + static void
94 + print_shuf(long *n, long aval)
95 + {
96 +        int     i;
97 +
98 +        for (i = 0; n[i+1]; i++) {
99 +                printf("\t%ld", aval % n[i]);
100 +                aval /= n[i];
101 +        }
102 +        printf("\t%ld\n", aval);
103 + }
104 +
105 +
106 + /* Allocate and prepare occupancy tree */
107 + static uby8 *
108 + tree_alloc(long alen)
109 + {
110 +        uby8  *troot;
111 +        double  bytes_per_bit;
112 +        int  i;
113 +        int  ht = 0;
114 +                                        /* how tall does our tree need to be? */
115 +        while (tree_br[ht].capacity*BRORDER < alen)
116 +                if (++ht >= NLEVELS) {
117 +                        fputs("Array too large to shuffle\n", stderr);
118 +                        exit(1);
119 +                }
120 +        bytes_per_bit = 1.;             /* figure out tree size (with overhead) */
121 +        for (i = ht; i >= 0; i--)
122 +                bytes_per_bit += (double)tree_br[i].cntr_siz;
123 +        bytes_per_bit += (double)tree_br[ht].skip_bytes;
124 +        bytes_per_bit /= (double)tree_br[ht].capacity;
125 +        troot = (uby8 *)calloc((long)(alen*bytes_per_bit)+2, 1);
126 +        if (troot == NULL) {
127 +                fputs("Not enough memory for shuffle\n", stderr);
128 +                exit(1);
129 +        }
130 +        *troot = ht;                    /* first byte is tree height */
131 +        for (i = 256; i--; ) {          /* assign 0-bit count table */
132 +                int     b;
133 +                buf[i] = 8;
134 +                for (b = i; b; b >>= 1)
135 +                        buf[i] -= b&1;
136 +        }
137 +        return(troot);
138 + }
139 +
140 +
141 + /* Get number of slots available at this branch location */
142 + static long
143 + get_avail(const uby8 *ctrp, int lvl)
144 + {
145 +        long    cnt = 0;
146 +        int     n = tree_br[lvl].cntr_siz;
147 +
148 +        while (--n > 0) {               /* LSB first */
149 +                cnt |= ctrp[n];
150 +                cnt <<= 8;
151 +        }
152 +        cnt |= ctrp[0];
153 +
154 +        return(tree_br[lvl].capacity - cnt);
155 + }
156 +
157 +
158 + /* Increment branch occupancy counter */
159 + static void
160 + incr_counter(uby8 *ctrp, int n)
161 + {
162 +        n = tree_br[n].cntr_siz;
163 +
164 +        while (++(*ctrp++))             /* LSB first */
165 +                if (--n <= 0) {
166 +                        fputs("Shuffle occupancy overflow!\n", stderr);
167 +                        exit(1);        /* means we sized something wrong */
168 +                }
169 + }
170 +
171 +
172 + /* Skip to and allocate a leaf from tree */
173 + static long
174 + eat_nth_leaf(uby8 *brp, long ski)
175 + {
176 +        int     lvl = *brp++;                   /* tree height in first byte */
177 +        long    pos = 0;
178 +        int     b;
179 +
180 +        while (lvl >= 0) {                      /* descend to leaves */
181 +                long  navail;
182 +                b = 0;                          /* select each branch */
183 +                while (ski >= (navail = get_avail(brp, lvl))) {
184 +                        if (++b >= BRORDER) {
185 +                                fputs("Shuffle tree error!\n", stderr);
186 +                                exit(1);
187 +                        }
188 +                        pos += tree_br[lvl].capacity;
189 +                        ski -= navail;
190 +                        brp += tree_br[lvl].skip_bytes;
191 +                }
192 +                incr_counter(brp, lvl);         /* we intend to eat one */
193 +                brp += tree_br[lvl--].cntr_siz; /* drop a level */
194 +        }
195 +        while (ski >= buf[*brp]) {              /* browse the leaves */
196 +                pos += 8;
197 +                ski -= buf[*brp++];             /* buf contains 0-bit counts */
198 +        }
199 +        b = 0;                                  /* find target bit in byte */
200 +        while ((ski -= !(*brp & 1<<b)) >= 0) {
201 +                pos++;
202 +                b++;
203 +        }
204 +        *brp |= 1<<b;                           /* eat it */
205 +        return(pos);                            /* & return leaf's slot# */
206 + }
207 +
208 +
209 + /* Shuffle all possible output strings and spit out randomly (tree version) */
210 + static void
211 + big_shuffle(long *n, long alen)
212 + {
213 +        uby8  *tree_root;
214 +                                        /* size and allocate holder tree */
215 +        tree_root = tree_alloc(alen);
216 +
217 +        while (alen > 0)                /* allocate and print random array entries */
218 +                print_shuf(n, eat_nth_leaf(tree_root, irandom(alen--)));
219 +
220 +        free(tree_root);                /* all done */
221 + }
222 +
223 +
224   /* Shuffle all possible output strings and spit out randomly */
225   static void
226 < shuffle(
74 < int *n
75 < )
226 > shuffle(long *n)
227   {
228 <        int  sub[MAXDIM];
229 <        int  ndim;
230 <        int  alen;
80 <        int  *myshuf;
81 <        int  i, j;
228 >        long  alen;
229 >        uint16  *myshuf;
230 >        int  i;
231  
232 <        alen = 1;               /* allocate shuffle index array */
233 <        for (j = 0; n[j]; j++)
234 <                if ((alen *= n[j]) < 0)
232 >        alen = 1;               /* compute shuffle size */
233 >        for (i = 0; n[i]; i++) {
234 >                if (alen*n[i] <= alen) {
235 >                        fputs("Array too large to count!\n", stderr);
236                          exit(1);
237 +                }
238 +                alen *= n[i];
239 +        }
240 +                                /* get unique starting point */
241 +        srandom((long)time(0));
242  
243 <        myshuf = (int *)malloc(alen*sizeof(int));
243 >        if (alen > 1L<<16) {    /* use large shuffle method? */
244 >                big_shuffle(n, alen);
245 >                return;
246 >        }
247 >        myshuf = (uint16 *)malloc(alen*sizeof(uint16));
248          if (myshuf == NULL) {
249                  fputs("Insufficient memory for shuffle\n", stderr);
250                  exit(1);
251          }
252          for (i = alen; i--; )   /* initialize in any order */
253                  myshuf[i] = i;
95 –                                /* get unique starting point */
96 –        srandom((long)time(0));
254                                  /* perform Fisher-Yates shuffle */
255          for (i = 0; i < alen-1; i++) {
256 <                int     ix = random()%(alen-i) + i;
256 >                int     ix = irandom(alen-i) + i;
257                  int     ndx = myshuf[i];
258                  myshuf[i] = myshuf[ix];
259                  myshuf[ix] = ndx;
260          }
261                                  /* put randomly indexed output */
262 <        for (i = alen; i--; ) {
263 <                int     aval = myshuf[i];
264 <                for (j = 0; n[j+1]; j++) {
265 <                        printf("\t%d", aval % n[j]);
109 <                        aval /= n[j];
110 <                }
111 <                printf("\t%d\n", aval);
112 <        }
113 <        free(myshuf);
262 >        for (i = alen; i--; )
263 >                print_shuf(n, (long)myshuf[i]);
264 >
265 >        free(myshuf);           /* all done */
266   }
267  
268  
269   int
270 < main(
119 < int  argc,
120 < char  *argv[]
121 < )
270 > main(int argc, char *argv[])
271   {
272          char  *prog = argv[0];
273          int  doshuffle = 0;
274 <        int  n[MAXDIM];
274 >        long  n[MAXDIM];
275          int  a;
276  
277          argv++; argc--;
# Line 133 | Line 282 | char  *argv[]
282                  argv++; argc--;
283          }
284          for (a = 0; a < argc; a++)
285 <                if ((n[a] = atoi(argv[a])) <= 1)
285 >                if ((n[a] = atol(argv[a])) <= 1)
286                          goto userr;
287          n[a] = 0;
288          if (!a)

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines