ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/hd/rholo3.c
Revision: 3.22
Committed: Wed Jul 8 17:59:58 1998 UTC (26 years, 3 months ago) by gwlarson
Content type: text/plain
Branch: MAIN
Changes since 3.21: +15 -18 lines
Log Message:
fixed bundle_set() with BS_ADJ to send satisfiable requests

File Contents

# Content
1 /*
2 * Routines for tracking beam compuatations
3 */
4
5 #include "rholo.h"
6
7 #define abs(x) ((x) > 0 ? (x) : -(x))
8 #define sgn(x) ((x) > 0 ? 1 : (x) < 0 ? -1 : 0)
9
10 static PACKHEAD *complist=NULL; /* list of beams to compute */
11 static int complen=0; /* length of complist */
12 static int listpos=0; /* current list position for next_packet */
13 static int lastin= -1; /* last ordered position in list */
14
15
16 int
17 beamcmp(b0, b1) /* comparison for descending compute order */
18 register PACKHEAD *b0, *b1;
19 {
20 return( b1->nr*(b0->nc+1) - b0->nr*(b1->nc+1) );
21 }
22
23
24 int
25 dispbeam(b, hp, bi) /* display a holodeck beam */
26 register BEAM *b;
27 HOLO *hp;
28 int bi;
29 {
30 static int n = 0;
31 static PACKHEAD *p = NULL;
32
33 if (b == NULL)
34 return;
35 if (b->nrm > n) { /* (re)allocate packet holder */
36 n = b->nrm;
37 if (p == NULL) p = (PACKHEAD *)malloc(packsiz(n));
38 else p = (PACKHEAD *)realloc((char *)p, packsiz(n));
39 if (p == NULL)
40 error(SYSTEM, "out of memory in dispbeam");
41 }
42 /* assign packet fields */
43 bcopy((char *)hdbray(b), (char *)packra(p), b->nrm*sizeof(RAYVAL));
44 p->nr = p->nc = b->nrm;
45 for (p->hd = 0; hdlist[p->hd] != hp; p->hd++)
46 if (hdlist[p->hd] == NULL)
47 error(CONSISTENCY, "unregistered holodeck in dispbeam");
48 p->bi = bi;
49 disp_packet(p); /* display it */
50 }
51
52
53 bundle_set(op, clist, nents) /* bundle set operation */
54 int op;
55 register PACKHEAD *clist;
56 int nents;
57 {
58 int oldnr;
59 register HDBEAMI *hb;
60 register int i, n;
61 /* look for common members */
62 for (n = 0; n < nents; n++) {
63 for (i = 0; i < complen; i++)
64 if (clist[n].bi == complist[i].bi &&
65 clist[n].hd == complist[i].hd) {
66 oldnr = complist[i].nr;
67 clist[n].nc = complist[i].nc;
68 switch (op) {
69 case BS_ADD: /* add to count */
70 complist[i].nr += clist[n].nr;
71 clist[n].nr = 0;
72 break;
73 case BS_ADJ: /* reset count */
74 complist[i].nr = clist[n].nr;
75 clist[n].nr = 0;
76 break;
77 case BS_DEL: /* delete count */
78 if (clist[n].nr == 0 ||
79 clist[n].nr >= complist[i].nr)
80 complist[i].nr = 0;
81 else
82 complist[i].nr -= clist[n].nr;
83 break;
84 }
85 if (complist[i].nr != oldnr)
86 lastin = -1; /* flag sort */
87 break;
88 }
89 if (i >= complen)
90 clist[n].nc = bnrays(hdlist[clist[n].hd], clist[n].bi);
91 }
92 /* complete list operations */
93 switch (op) {
94 case BS_NEW: /* new computation set */
95 listpos = 0; lastin = -1;
96 if (complen) /* free old list */
97 free((char *)complist);
98 complist = NULL;
99 if (!(complen = nents))
100 return;
101 complist = (PACKHEAD *)malloc(nents*sizeof(PACKHEAD));
102 if (complist == NULL)
103 goto memerr;
104 bcopy((char *)clist, (char *)complist, nents*sizeof(PACKHEAD));
105 break;
106 case BS_ADD: /* add to computation set */
107 case BS_ADJ: /* adjust set quantities */
108 if (nents <= 0)
109 return;
110 sortcomplist(); /* sort updated list & new entries */
111 qsort((char *)clist, nents, sizeof(PACKHEAD), beamcmp);
112 /* what can't we satisfy? */
113 for (n = 0; n < nents && clist[n].nr > clist[n].nc; n++)
114 ;
115 if (op == BS_ADJ) { /* don't regenerate adjusted beams */
116 for (i = n; i < nents && clist[i].nr > 0; i++)
117 ;
118 nents = i;
119 }
120 if (n) { /* allocate space for merged list */
121 PACKHEAD *newlist;
122 newlist = (PACKHEAD *)malloc(
123 (complen+n)*sizeof(PACKHEAD) );
124 if (newlist == NULL)
125 goto memerr;
126 /* merge lists */
127 mergeclists(newlist, clist, n, complist, complen);
128 if (complen)
129 free((char *)complist);
130 complist = newlist;
131 complen += n;
132 }
133 listpos = 0;
134 lastin = complen-1; /* list is now sorted */
135 break;
136 case BS_DEL: /* delete from computation set */
137 return; /* already done */
138 default:
139 error(CONSISTENCY, "bundle_set called with unknown operation");
140 }
141 if (outdev == NULL) /* nothing to display? */
142 return;
143 /* load and display beams we have */
144 hb = (HDBEAMI *)malloc(nents*sizeof(HDBEAMI));
145 for (i = 0; i < nents; i++) {
146 hb[i].h = hdlist[clist[i].hd];
147 hb[i].b = clist[i].bi;
148 }
149 hdloadbeams(hb, nents, dispbeam);
150 free((char *)hb);
151 return;
152 memerr:
153 error(SYSTEM, "out of memory in bundle_set");
154 }
155
156
157 double
158 beamvolume(hp, bi) /* compute approximate volume of a beam */
159 HOLO *hp;
160 int bi;
161 {
162 GCOORD gc[2];
163 FVECT cp[4], edgeA, edgeB, cent[2];
164 FVECT v, crossp[2], diffv;
165 double vol[2];
166 register int i;
167 /* get grid coordinates */
168 if (!hdbcoord(gc, hp, bi))
169 error(CONSISTENCY, "bad beam index in beamvolume");
170 for (i = 0; i < 2; i++) { /* compute cell area vectors */
171 hdcell(cp, hp, gc+i);
172 VSUM(edgeA, cp[1], cp[0], -1.0);
173 VSUM(edgeB, cp[3], cp[1], -1.0);
174 fcross(crossp[i], edgeA, edgeB);
175 /* compute center */
176 cent[i][0] = 0.5*(cp[0][0] + cp[2][0]);
177 cent[i][1] = 0.5*(cp[0][1] + cp[2][1]);
178 cent[i][2] = 0.5*(cp[0][2] + cp[2][2]);
179 }
180 /* compute difference vector */
181 VSUM(diffv, cent[1], cent[0], -1.0);
182 for (i = 0; i < 2; i++) { /* compute volume contributions */
183 vol[i] = 0.5*DOT(crossp[i], diffv);
184 if (vol[i] < 0.) vol[i] = -vol[i];
185 }
186 return(vol[0] + vol[1]); /* return total volume */
187 }
188
189
190 init_global() /* initialize global ray computation */
191 {
192 long wtotal = 0;
193 double frac;
194 int i;
195 register int j, k;
196 /* free old list and empty queue */
197 if (complen > 0) {
198 free((char *)complist);
199 done_packets(flush_queue());
200 }
201 /* allocate beam list */
202 complen = 0;
203 for (j = 0; hdlist[j] != NULL; j++)
204 complen += nbeams(hdlist[j]);
205 complist = (PACKHEAD *)malloc(complen*sizeof(PACKHEAD));
206 if (complist == NULL)
207 error(SYSTEM, "out of memory in init_global");
208 /* compute beam weights */
209 k = 0;
210 for (j = 0; hdlist[j] != NULL; j++) {
211 frac = 512. * VLEN(hdlist[j]->wg[0]) *
212 VLEN(hdlist[j]->wg[1]) *
213 VLEN(hdlist[j]->wg[2]);
214 for (i = nbeams(hdlist[j]); i > 0; i--) {
215 complist[k].hd = j;
216 complist[k].bi = i;
217 complist[k].nr = frac*beamvolume(hdlist[j], i) + 0.5;
218 complist[k].nc = bnrays(hdlist[j], i);
219 wtotal += complist[k++].nr;
220 }
221 }
222 /* adjust weights */
223 if (vdef(DISKSPACE))
224 frac = 1024.*1024.*vflt(DISKSPACE) / (wtotal*sizeof(RAYVAL));
225 else
226 frac = 1024.*1024.*16384. / (wtotal*sizeof(RAYVAL));
227 while (k--)
228 complist[k].nr = frac * complist[k].nr;
229 listpos = 0; lastin = -1; /* perform initial sort */
230 sortcomplist();
231 }
232
233
234 mergeclists(cdest, cl1, n1, cl2, n2) /* merge two sorted lists */
235 register PACKHEAD *cdest;
236 register PACKHEAD *cl1, *cl2;
237 int n1, n2;
238 {
239 register int cmp;
240
241 while (n1 | n2) {
242 if (!n1) cmp = 1;
243 else if (!n2) cmp = -1;
244 else cmp = beamcmp(cl1, cl2);
245 if (cmp > 0) {
246 copystruct(cdest, cl2);
247 cl2++; n2--;
248 } else {
249 copystruct(cdest, cl1);
250 cl1++; n1--;
251 }
252 cdest++;
253 }
254 }
255
256
257 sortcomplist() /* fix our list order */
258 {
259 PACKHEAD *list2;
260 register int i;
261
262 if (complen <= 0) /* check to see if there is even a list */
263 return;
264 if (lastin < 0 || listpos*4 >= complen*3)
265 qsort((char *)complist, complen, sizeof(PACKHEAD), beamcmp);
266 else if (listpos) { /* else sort and merge sublist */
267 list2 = (PACKHEAD *)malloc(listpos*sizeof(PACKHEAD));
268 if (list2 == NULL)
269 error(SYSTEM, "out of memory in sortcomplist");
270 bcopy((char *)complist,(char *)list2,listpos*sizeof(PACKHEAD));
271 qsort((char *)list2, listpos, sizeof(PACKHEAD), beamcmp);
272 mergeclists(complist, list2, listpos,
273 complist+listpos, complen-listpos);
274 free((char *)list2);
275 }
276 /* drop satisfied requests */
277 for (i = complen; i-- && complist[i].nr <= complist[i].nc; )
278 ;
279 if (i < 0) {
280 free((char *)complist);
281 complist = NULL;
282 complen = 0;
283 } else if (i < complen-1) {
284 list2 = (PACKHEAD *)realloc((char *)complist,
285 (i+1)*sizeof(PACKHEAD));
286 if (list2 != NULL) {
287 complist = list2;
288 complen = i+1;
289 }
290 }
291 listpos = 0; lastin = i;
292 }
293
294
295 /*
296 * The following routine works on the assumption that the bundle weights are
297 * more or less evenly distributed, such that computing a packet causes
298 * a given bundle to move way down in the computation order. We keep
299 * track of where the computed bundle with the highest priority would end
300 * up, and if we get further in our compute list than this, we resort the
301 * list and start again from the beginning. Since
302 * a merge sort is used, the sorting costs are minimal.
303 */
304 next_packet(p, n) /* prepare packet for computation */
305 register PACKET *p;
306 int n;
307 {
308 register int i;
309
310 if (listpos > lastin) /* time to sort the list */
311 sortcomplist();
312 if (complen <= 0)
313 return(0);
314 p->hd = complist[listpos].hd;
315 p->bi = complist[listpos].bi;
316 p->nc = complist[listpos].nc;
317 p->nr = complist[listpos].nr - p->nc;
318 if (p->nr <= 0)
319 return(0);
320 #ifdef DEBUG
321 if (n < 1 | n > RPACKSIZ)
322 error(CONSISTENCY, "next_packet called with bad n value");
323 #endif
324 if (p->nr > n)
325 p->nr = n;
326 complist[listpos].nc += p->nr; /* find where this one would go */
327 while (lastin > listpos &&
328 beamcmp(complist+lastin, complist+listpos) > 0)
329 lastin--;
330 listpos++;
331 return(1);
332 }