| 10 |  |  | 
| 11 |  | #include "rholo.h" | 
| 12 |  |  | 
| 13 | – |  | 
| 13 |  | #define abs(x)          ((x) > 0 ? (x) : -(x)) | 
| 14 |  | #define sgn(x)          ((x) > 0 ? 1 : (x) < 0 ? -1 : 0) | 
| 15 |  |  | 
| 17 | – |  | 
| 16 |  | static PACKHEAD *complist=NULL; /* list of beams to compute */ | 
| 17 |  | static int      complen=0;      /* length of complist */ | 
| 18 |  | static int      listpos=0;      /* current list position for next_packet */ | 
| 23 |  | beamcmp(b0, b1)                 /* comparison for descending compute order */ | 
| 24 |  | register PACKHEAD       *b0, *b1; | 
| 25 |  | { | 
| 26 | < | return( b1->nr*(bnrays(hdlist[b0->hd],b0->bi)+1) - | 
| 29 | < | b0->nr*(bnrays(hdlist[b1->hd],b1->bi)+1) ); | 
| 26 | > | return( b1->nr*(b0->nc+1) - b0->nr*(b1->nc+1) ); | 
| 27 |  | } | 
| 28 |  |  | 
| 29 |  |  | 
| 30 | + | int | 
| 31 | + | dispbeam(b, hp, bi)                     /* display a holodeck beam */ | 
| 32 | + | register BEAM   *b; | 
| 33 | + | HOLO    *hp; | 
| 34 | + | int     bi; | 
| 35 | + | { | 
| 36 | + | static int      n = 0; | 
| 37 | + | static PACKHEAD *p = NULL; | 
| 38 | + |  | 
| 39 | + | if (b == NULL) | 
| 40 | + | return; | 
| 41 | + | if (b->nrm > n) {               /* (re)allocate packet holder */ | 
| 42 | + | n = b->nrm; | 
| 43 | + | if (p == NULL) p = (PACKHEAD *)malloc(packsiz(n)); | 
| 44 | + | else p = (PACKHEAD *)realloc((char *)p, packsiz(n)); | 
| 45 | + | if (p == NULL) | 
| 46 | + | error(SYSTEM, "out of memory in dispbeam"); | 
| 47 | + | } | 
| 48 | + | /* assign packet fields */ | 
| 49 | + | bcopy((char *)hdbray(b), (char *)packra(p), b->nrm*sizeof(RAYVAL)); | 
| 50 | + | p->nr = p->nc = b->nrm; | 
| 51 | + | for (p->hd = 0; hdlist[p->hd] != hp; p->hd++) | 
| 52 | + | if (hdlist[p->hd] == NULL) | 
| 53 | + | error(CONSISTENCY, "unregistered holodeck in dispbeam"); | 
| 54 | + | p->bi = bi; | 
| 55 | + | disp_packet(p);                 /* display it */ | 
| 56 | + | } | 
| 57 | + |  | 
| 58 | + |  | 
| 59 |  | bundle_set(op, clist, nents)    /* bundle set operation */ | 
| 60 |  | int     op; | 
| 61 | < | PACKHEAD        *clist; | 
| 61 | > | register PACKHEAD       *clist; | 
| 62 |  | int     nents; | 
| 63 |  | { | 
| 38 | – | BEAM    *b; | 
| 39 | – | PACKHEAD        *p; | 
| 64 |  | register int    i, n; | 
| 65 |  |  | 
| 66 |  | switch (op) { | 
| 77 |  | if (complist == NULL) | 
| 78 |  | goto memerr; | 
| 79 |  | bcopy((char *)clist, (char *)complist, nents*sizeof(PACKHEAD)); | 
| 80 | < | complen = nents; | 
| 57 | < | listpos = 0; | 
| 58 | < | lastin = -1;            /* flag for initial sort */ | 
| 80 | > | complen = nents;        /* finish initialization below */ | 
| 81 |  | break; | 
| 82 |  | case BS_ADD:                    /* add to computation set */ | 
| 83 | + | case BS_ADJ:                    /* adjust set quantities */ | 
| 84 |  | if (nents <= 0) | 
| 85 |  | return; | 
| 86 |  | /* merge any common members */ | 
| 87 | < | for (i = 0; i < complen; i++) | 
| 88 | < | for (n = 0; n < nents; n++) | 
| 87 | > | for (n = 0; n < nents; n++) { | 
| 88 | > | for (i = 0; i < complen; i++) | 
| 89 |  | if (clist[n].bi == complist[i].bi && | 
| 90 |  | clist[n].hd == complist[i].hd) { | 
| 91 | < | complist[i].nr += clist[n].nr; | 
| 91 | > | int     oldnr = complist[i].nr; | 
| 92 | > | if (op == BS_ADD) | 
| 93 | > | complist[i].nr += clist[n].nr; | 
| 94 | > | else /* op == BS_ADJ */ | 
| 95 | > | complist[i].nr = clist[n].nr; | 
| 96 |  | clist[n].nr = 0; | 
| 97 | < | lastin = -1;    /* flag full sort */ | 
| 97 | > | clist[n].nc = complist[i].nc; | 
| 98 | > | if (complist[i].nr != oldnr) | 
| 99 | > | lastin = -1;    /* flag sort */ | 
| 100 |  | break; | 
| 101 |  | } | 
| 102 | + | if (i >= complen) | 
| 103 | + | clist[n].nc = bnrays(hdlist[clist[n].hd], | 
| 104 | + | clist[n].bi); | 
| 105 | + | } | 
| 106 |  | /* sort updated list */ | 
| 107 |  | sortcomplist(); | 
| 108 |  | /* sort new entries */ | 
| 109 |  | qsort((char *)clist, nents, sizeof(PACKHEAD), beamcmp); | 
| 110 |  | /* what can't we satisfy? */ | 
| 111 | < | for (n = 0; n < nents && clist[n].nr > | 
| 79 | < | bnrays(hdlist[clist[n].hd],clist[n].bi); n++) | 
| 111 | > | for (n = 0; n < nents && clist[n].nr > clist[n].nc; n++) | 
| 112 |  | ; | 
| 113 | + | if (op == BS_ADJ) | 
| 114 | + | nents = n; | 
| 115 |  | if (n) {                /* allocate space for merged list */ | 
| 116 |  | PACKHEAD        *newlist; | 
| 117 |  | newlist = (PACKHEAD *)malloc( | 
| 144 |  | lastin = -1;    /* flag full sort */ | 
| 145 |  | break; | 
| 146 |  | } | 
| 113 | – | if (lastin < 0)         /* sort updated list */ | 
| 114 | – | sortcomplist(); | 
| 147 |  | return;                 /* no display */ | 
| 148 |  | default: | 
| 149 |  | error(CONSISTENCY, "bundle_set called with unknown operation"); | 
| 150 |  | } | 
| 151 | < | n = 0;                          /* allocate packet holder */ | 
| 152 | < | for (i = 0; i < nents; i++) | 
| 153 | < | if (clist[i].nr > n) | 
| 154 | < | n = clist[i].nr; | 
| 155 | < | p = (PACKHEAD *)malloc(sizeof(PACKHEAD) + n*sizeof(RAYVAL)); | 
| 156 | < | if (p == NULL) | 
| 157 | < | goto memerr; | 
| 126 | < | /* display what we have */ | 
| 127 | < | for (i = 0; i < nents; i++) | 
| 128 | < | if ((b = hdgetbeam(hdlist[clist[i].hd], clist[i].bi)) != NULL) { | 
| 129 | < | bcopy((char *)hdbray(b), (char *)(p+1), | 
| 130 | < | (p->nr=b->nrm)*sizeof(RAYVAL)); | 
| 131 | < | disp_packet((PACKET *)p); | 
| 151 | > | if (outdev != NULL) {           /* load and display beams we have */ | 
| 152 | > | register HDBEAMI        *hb; | 
| 153 | > |  | 
| 154 | > | hb = (HDBEAMI *)malloc(nents*sizeof(HDBEAMI)); | 
| 155 | > | for (i = 0; i < nents; i++) { | 
| 156 | > | hb[i].h = hdlist[clist[i].hd]; | 
| 157 | > | hb[i].b = clist[i].bi; | 
| 158 |  | } | 
| 159 | < | free((char *)p);                /* clean up */ | 
| 159 | > | hdloadbeams(hb, nents, dispbeam); | 
| 160 | > | free((char *)hb); | 
| 161 | > | } | 
| 162 | > | if (op == BS_NEW) { | 
| 163 | > | done_packets(flush_queue());    /* empty queue, so we can... */ | 
| 164 | > | for (i = 0; i < complen; i++)   /* ...get number computed */ | 
| 165 | > | complist[i].nc = bnrays(hdlist[complist[i].hd], | 
| 166 | > | complist[i].bi); | 
| 167 | > | listpos = 0; | 
| 168 | > | lastin = -1;            /* flag for initial sort */ | 
| 169 | > | } | 
| 170 |  | return; | 
| 171 |  | memerr: | 
| 172 |  | error(SYSTEM, "out of memory in bundle_set"); | 
| 173 |  | } | 
| 174 |  |  | 
| 175 |  |  | 
| 176 | < | int | 
| 177 | < | weightf(hp, x0, x1, x2)         /* voxel weighting function */ | 
| 142 | < | register HOLO   *hp; | 
| 143 | < | register int    x0, x1, x2; | 
| 144 | < | { | 
| 145 | < | switch (vlet(OCCUPANCY)) { | 
| 146 | < | case 'U':               /* uniform weighting */ | 
| 147 | < | return(1); | 
| 148 | < | case 'C':               /* center weighting (crude) */ | 
| 149 | < | x0 += x0 - hp->grid[0] + 1; | 
| 150 | < | x0 = abs(x0)*hp->grid[1]*hp->grid[2]; | 
| 151 | < | x1 += x1 - hp->grid[1] + 1; | 
| 152 | < | x1 = abs(x1)*hp->grid[0]*hp->grid[2]; | 
| 153 | < | x2 += x2 - hp->grid[2] + 1; | 
| 154 | < | x2 = abs(x2)*hp->grid[0]*hp->grid[1]; | 
| 155 | < | return(hp->grid[0]*hp->grid[1]*hp->grid[2] - | 
| 156 | < | (x0+x1+x2)/3); | 
| 157 | < | default: | 
| 158 | < | badvalue(OCCUPANCY); | 
| 159 | < | } | 
| 160 | < | } | 
| 161 | < |  | 
| 162 | < |  | 
| 163 | < | /* The following is by Daniel Cohen, taken from Graphics Gems IV, p. 368 */ | 
| 164 | < | long | 
| 165 | < | lineweight(hp, x, y, z, dx, dy, dz)     /* compute weights along a line */ | 
| 176 | > | double | 
| 177 | > | beamvolume(hp, bi)      /* compute approximate volume of a beam */ | 
| 178 |  | HOLO    *hp; | 
| 179 | < | int     x, y, z, dx, dy, dz; | 
| 179 | > | int     bi; | 
| 180 |  | { | 
| 181 | < | long    wres = 0; | 
| 182 | < | int     n, sx, sy, sz, exy, exz, ezy, ax, ay, az, bx, by, bz; | 
| 183 | < |  | 
| 184 | < | sx = sgn(dx);   sy = sgn(dy);   sz = sgn(dz); | 
| 185 | < | ax = abs(dx);   ay = abs(dy);   az = abs(dz); | 
| 186 | < | bx = 2*ax;      by = 2*ay;      bz = 2*az; | 
| 187 | < | exy = ay-ax;    exz = az-ax;    ezy = ay-az; | 
| 188 | < | n = ax+ay+az + 1;               /* added increment to visit last */ | 
| 189 | < | while (n--) { | 
| 190 | < | wres += weightf(hp, x, y, z); | 
| 191 | < | if (exy < 0) { | 
| 192 | < | if (exz < 0) { | 
| 193 | < | x += sx; | 
| 194 | < | exy += by; exz += bz; | 
| 195 | < | } else { | 
| 196 | < | z += sz; | 
| 197 | < | exz -= bx; ezy += by; | 
| 186 | < | } | 
| 187 | < | } else { | 
| 188 | < | if (ezy < 0) { | 
| 189 | < | z += sz; | 
| 190 | < | exz -= bx; ezy += by; | 
| 191 | < | } else { | 
| 192 | < | y += sy; | 
| 193 | < | exy -= bx; ezy -= bz; | 
| 194 | < | } | 
| 195 | < | } | 
| 181 | > | GCOORD  gc[2]; | 
| 182 | > | FVECT   cp[4], edgeA, edgeB, cent[2]; | 
| 183 | > | FVECT   v, crossp[2], diffv; | 
| 184 | > | double  vol[2]; | 
| 185 | > | register int    i; | 
| 186 | > | /* get grid coordinates */ | 
| 187 | > | if (!hdbcoord(gc, hp, bi)) | 
| 188 | > | error(CONSISTENCY, "bad beam index in beamvolume"); | 
| 189 | > | for (i = 0; i < 2; i++) {       /* compute cell area vectors */ | 
| 190 | > | hdcell(cp, hp, gc+i); | 
| 191 | > | VSUM(edgeA, cp[1], cp[0], -1.0); | 
| 192 | > | VSUM(edgeB, cp[3], cp[1], -1.0); | 
| 193 | > | fcross(crossp[i], edgeA, edgeB); | 
| 194 | > | /* compute center */ | 
| 195 | > | cent[i][0] = 0.5*(cp[0][0] + cp[2][0]); | 
| 196 | > | cent[i][1] = 0.5*(cp[0][1] + cp[2][1]); | 
| 197 | > | cent[i][2] = 0.5*(cp[0][2] + cp[2][2]); | 
| 198 |  | } | 
| 199 | < | return(wres); | 
| 199 | > | /* compute difference vector */ | 
| 200 | > | VSUM(diffv, cent[1], cent[0], -1.0); | 
| 201 | > | for (i = 0; i < 2; i++) {       /* compute volume contributions */ | 
| 202 | > | vol[i] = 0.5*DOT(crossp[i], diffv); | 
| 203 | > | if (vol[i] < 0.) vol[i] = -vol[i]; | 
| 204 | > | } | 
| 205 | > | return(vol[0] + vol[1]);        /* return total volume */ | 
| 206 |  | } | 
| 207 |  |  | 
| 208 |  |  | 
| 209 |  | init_global()                   /* initialize global ray computation */ | 
| 210 |  | { | 
| 211 |  | long    wtotal = 0; | 
| 204 | – | int     i, j; | 
| 205 | – | int     lseg[2][3]; | 
| 212 |  | double  frac; | 
| 213 | < | register int    k; | 
| 213 | > | int     i; | 
| 214 | > | register int    j, k; | 
| 215 |  | /* free old list */ | 
| 216 |  | if (complen > 0) | 
| 217 |  | free((char *)complist); | 
| 224 |  | error(SYSTEM, "out of memory in init_global"); | 
| 225 |  | /* compute beam weights */ | 
| 226 |  | k = 0; | 
| 227 | < | for (j = 0; hdlist[j] != NULL; j++) | 
| 227 | > | for (j = 0; hdlist[j] != NULL; j++) { | 
| 228 | > | frac = 512. * hdlist[j]->wg[0] * | 
| 229 | > | hdlist[j]->wg[1] * hdlist[j]->wg[2]; | 
| 230 | > | if (frac < 0.) frac = -frac; | 
| 231 |  | for (i = nbeams(hdlist[j]); i > 0; i--) { | 
| 222 | – | hdlseg(lseg, hdlist[j], i); | 
| 232 |  | complist[k].hd = j; | 
| 233 |  | complist[k].bi = i; | 
| 234 | < | complist[k].nr = lineweight( hdlist[j], | 
| 226 | < | lseg[0][0], lseg[0][1], lseg[0][2], | 
| 227 | < | lseg[1][0] - lseg[0][0], | 
| 228 | < | lseg[1][1] - lseg[0][1], | 
| 229 | < | lseg[1][2] - lseg[0][2] ); | 
| 234 | > | complist[k].nr = frac*beamvolume(hdlist[j], i) + 0.5; | 
| 235 |  | wtotal += complist[k++].nr; | 
| 236 |  | } | 
| 237 | + | } | 
| 238 |  | /* adjust weights */ | 
| 239 | < | if (vdef(DISKSPACE)) { | 
| 239 | > | if (vdef(DISKSPACE)) | 
| 240 |  | frac = 1024.*1024.*vflt(DISKSPACE) / (wtotal*sizeof(RAYVAL)); | 
| 241 | < | if (frac < 0.95 | frac > 1.05) | 
| 242 | < | while (k--) | 
| 243 | < | complist[k].nr = frac * complist[k].nr; | 
| 244 | < | } | 
| 241 | > | else | 
| 242 | > | frac = 1024.*1024.*16384. / (wtotal*sizeof(RAYVAL)); | 
| 243 | > | while (k--) | 
| 244 | > | complist[k].nr = frac * complist[k].nr; | 
| 245 |  | listpos = 0; lastin = -1;       /* flag initial sort */ | 
| 246 |  | } | 
| 247 |  |  | 
| 248 |  |  | 
| 249 |  | mergeclists(cdest, cl1, n1, cl2, n2)    /* merge two sorted lists */ | 
| 250 | < | PACKHEAD        *cdest; | 
| 251 | < | PACKHEAD        *cl1, *cl2; | 
| 250 | > | register PACKHEAD       *cdest; | 
| 251 | > | register PACKHEAD       *cl1, *cl2; | 
| 252 |  | int     n1, n2; | 
| 253 |  | { | 
| 254 | < | int     cmp; | 
| 254 | > | register int    cmp; | 
| 255 |  |  | 
| 256 |  | while (n1 | n2) { | 
| 257 |  | if (!n1) cmp = 1; | 
| 274 |  | PACKHEAD        *list2; | 
| 275 |  | register int    i; | 
| 276 |  |  | 
| 271 | – | /* empty queue */ | 
| 272 | – | done_packets(flush_queue()); | 
| 277 |  | if (complen <= 0)       /* check to see if there is even a list */ | 
| 278 |  | return; | 
| 279 | < | if (lastin < 0)         /* flag to sort entire list */ | 
| 279 | > | if (lastin < 0 || listpos*4 >= complen*3) | 
| 280 |  | qsort((char *)complist, complen, sizeof(PACKHEAD), beamcmp); | 
| 281 |  | else if (listpos) {     /* else sort and merge sublist */ | 
| 282 |  | list2 = (PACKHEAD *)malloc(listpos*sizeof(PACKHEAD)); | 
| 288 |  | complist+listpos, complen-listpos); | 
| 289 |  | free((char *)list2); | 
| 290 |  | } | 
| 287 | – | /* check for all finished */ | 
| 288 | – | if (complist[0].nr <= bnrays(hdlist[complist[0].hd],complist[0].bi)) { | 
| 289 | – | free((char *)complist); | 
| 290 | – | complist = NULL; | 
| 291 | – | complen = 0; | 
| 292 | – | } | 
| 291 |  | /* drop satisfied requests */ | 
| 292 | < | for (i = complen; i-- && complist[i].nr <= | 
| 295 | < | bnrays(hdlist[complist[i].hd],complist[i].bi); ) | 
| 292 | > | for (i = complen; i-- && complist[i].nr <= complist[i].nc; ) | 
| 293 |  | ; | 
| 294 |  | if (i < 0) { | 
| 295 |  | free((char *)complist); | 
| 313 |  | * a given bundle to move way down in the computation order.  We keep | 
| 314 |  | * track of where the computed bundle with the highest priority would end | 
| 315 |  | * up, and if we get further in our compute list than this, we resort the | 
| 316 | < | * list and start again from the beginning.  We have to flush the queue | 
| 317 | < | * each time we sort, to ensure that we are not disturbing the order. | 
| 321 | < | *      If our major assumption is violated, and we have a very steep | 
| 322 | < | * descent in our weights, then we will end up resorting much more often | 
| 323 | < | * than necessary, resulting in frequent flushing of the queue.  Since | 
| 324 | < | * a merge sort is used, the sorting costs will be minimal. | 
| 316 | > | * list and start again from the beginning.  Since | 
| 317 | > | * a merge sort is used, the sorting costs are minimal. | 
| 318 |  | */ | 
| 319 | < | next_packet(p)                  /* prepare packet for computation */ | 
| 319 | > | next_packet(p, n)               /* prepare packet for computation */ | 
| 320 |  | register PACKET *p; | 
| 321 | + | int     n; | 
| 322 |  | { | 
| 329 | – | int     ncomp; | 
| 323 |  | register int    i; | 
| 324 |  |  | 
| 332 | – | if (complen <= 0) | 
| 333 | – | return(0); | 
| 325 |  | if (listpos > lastin)           /* time to sort the list */ | 
| 326 |  | sortcomplist(); | 
| 327 | + | if (complen <= 0) | 
| 328 | + | return(0); | 
| 329 |  | p->hd = complist[listpos].hd; | 
| 330 |  | p->bi = complist[listpos].bi; | 
| 331 | < | ncomp = bnrays(hdlist[p->hd],p->bi); | 
| 332 | < | p->nr = complist[listpos].nr - ncomp; | 
| 331 | > | p->nc = complist[listpos].nc; | 
| 332 | > | p->nr = complist[listpos].nr - p->nc; | 
| 333 |  | if (p->nr <= 0) | 
| 334 |  | return(0); | 
| 335 | < | if (p->nr > RPACKSIZ) | 
| 336 | < | p->nr = RPACKSIZ; | 
| 337 | < | ncomp += p->nr;                 /* find where this one would go */ | 
| 338 | < | while (lastin > listpos && complist[listpos].nr * | 
| 339 | < | (bnrays(hdlist[complist[lastin].hd],complist[lastin].bi)+1) | 
| 340 | < | > complist[lastin].nr * (ncomp+1)) | 
| 335 | > | #ifdef DEBUG | 
| 336 | > | if (n < 1 | n > RPACKSIZ) | 
| 337 | > | error(CONSISTENCY, "next_packet called with bad n value"); | 
| 338 | > | #endif | 
| 339 | > | if (p->nr > n) | 
| 340 | > | p->nr = n; | 
| 341 | > | complist[listpos].nc += p->nr;  /* find where this one would go */ | 
| 342 | > | while (lastin > listpos && | 
| 343 | > | beamcmp(complist+lastin, complist+listpos) > 0) | 
| 344 |  | lastin--; | 
| 345 |  | listpos++; | 
| 346 |  | return(1); |