| 1 | 
– | 
/* Copyright (c) 1998 Silicon Graphics, Inc. */ | 
| 2 | 
– | 
 | 
| 1 | 
  | 
#ifndef lint | 
| 2 | 
< | 
static char SCCSid[] = "$SunId$ SGI"; | 
| 2 | 
> | 
static const char       RCSid[] = "$Id$"; | 
| 3 | 
  | 
#endif | 
| 6 | 
– | 
 | 
| 4 | 
  | 
/* | 
| 5 | 
  | 
 * Rtrace support routines for holodeck rendering | 
| 6 | 
  | 
 */ | 
| 12 | 
  | 
 | 
| 13 | 
  | 
VIEWPOINT       myeye;          /* target view position */ | 
| 14 | 
  | 
 | 
| 15 | 
+ | 
struct gclim { | 
| 16 | 
+ | 
        HOLO    *hp;                    /* holodeck pointer */ | 
| 17 | 
+ | 
        GCOORD  gc;                     /* grid cell */ | 
| 18 | 
+ | 
        FVECT   egp;                    /* eye grid point */ | 
| 19 | 
+ | 
        double  erg2;                   /* mean square eye grid range */ | 
| 20 | 
+ | 
        double  gmin[2], gmax[2];       /* grid coordinate limits */ | 
| 21 | 
+ | 
};                              /* a grid coordinate range */ | 
| 22 | 
  | 
 | 
| 23 | 
+ | 
 | 
| 24 | 
+ | 
static | 
| 25 | 
+ | 
initeyelim(gcl, hp, gc)         /* initialize grid coordinate limits */ | 
| 26 | 
+ | 
register struct gclim   *gcl; | 
| 27 | 
+ | 
register HOLO   *hp; | 
| 28 | 
+ | 
GCOORD  *gc; | 
| 29 | 
+ | 
{ | 
| 30 | 
+ | 
        register RREAL  *v; | 
| 31 | 
+ | 
        register int    i; | 
| 32 | 
+ | 
 | 
| 33 | 
+ | 
        if (hp != NULL) { | 
| 34 | 
+ | 
                hdgrid(gcl->egp, gcl->hp = hp, myeye.vpt); | 
| 35 | 
+ | 
                gcl->erg2 = 0; | 
| 36 | 
+ | 
                for (i = 0, v = hp->wg[0]; i < 3; i++, v += 3) | 
| 37 | 
+ | 
                        gcl->erg2 += DOT(v,v); | 
| 38 | 
+ | 
                gcl->erg2 *= (1./3.) * myeye.rng*myeye.rng; | 
| 39 | 
+ | 
        } | 
| 40 | 
+ | 
        if (gc != NULL) | 
| 41 | 
+ | 
                copystruct(&gcl->gc, gc); | 
| 42 | 
+ | 
        gcl->gmin[0] = gcl->gmin[1] = FHUGE; | 
| 43 | 
+ | 
        gcl->gmax[0] = gcl->gmax[1] = -FHUGE; | 
| 44 | 
+ | 
} | 
| 45 | 
+ | 
 | 
| 46 | 
+ | 
 | 
| 47 | 
+ | 
static | 
| 48 | 
+ | 
groweyelim(gcl, gc, r0, r1, tight)      /* grow grid limits about eye point */ | 
| 49 | 
+ | 
register struct gclim   *gcl; | 
| 50 | 
+ | 
GCOORD  *gc; | 
| 51 | 
+ | 
double  r0, r1; | 
| 52 | 
+ | 
int     tight; | 
| 53 | 
+ | 
{ | 
| 54 | 
+ | 
        FVECT   gp, ab; | 
| 55 | 
+ | 
        double  ab2, od, cfact; | 
| 56 | 
+ | 
        double  sqcoef[3], ctcoef[3], licoef[3], cnst; | 
| 57 | 
+ | 
        int     gw, gi[2]; | 
| 58 | 
+ | 
        double  wallpos, a, b, c, d, e, f; | 
| 59 | 
+ | 
        double  root[2], yex; | 
| 60 | 
+ | 
        int     n, i, j, nex; | 
| 61 | 
+ | 
                                                /* point/view cone */ | 
| 62 | 
+ | 
        i = gc->w>>1; | 
| 63 | 
+ | 
        gp[i] = gc->w&1 ? gcl->hp->grid[i] : 0; | 
| 64 | 
+ | 
        gp[hdwg0[gc->w]] = gc->i[0] + r0; | 
| 65 | 
+ | 
        gp[hdwg1[gc->w]] = gc->i[1] + r1; | 
| 66 | 
+ | 
        VSUB(ab, gcl->egp, gp); | 
| 67 | 
+ | 
        ab2 = DOT(ab, ab); | 
| 68 | 
+ | 
        gw = gcl->gc.w>>1; | 
| 69 | 
+ | 
        if ((i==gw ? ab[gw]*ab[gw] : ab2)  <= gcl->erg2 + FTINY) { | 
| 70 | 
+ | 
                gcl->gmin[0] = gcl->gmin[1] = -FHUGE; | 
| 71 | 
+ | 
                gcl->gmax[0] = gcl->gmax[1] = FHUGE; | 
| 72 | 
+ | 
                return;                 /* too close (to wall) */ | 
| 73 | 
+ | 
        } | 
| 74 | 
+ | 
        ab2 = 1./ab2;                           /* 1/norm2(ab) */ | 
| 75 | 
+ | 
        od = DOT(gp, ab);                       /* origin dot direction */ | 
| 76 | 
+ | 
        cfact = 1./(1. - ab2*gcl->erg2);        /* tan^2 + 1 of cone angle */ | 
| 77 | 
+ | 
        for (i = 0; i < 3; i++) {               /* compute cone equation */ | 
| 78 | 
+ | 
                sqcoef[i] = ab[i]*ab[i]*cfact*ab2 - 1.; | 
| 79 | 
+ | 
                ctcoef[i] = 2.*ab[i]*ab[(i+1)%3]*cfact*ab2; | 
| 80 | 
+ | 
                licoef[i] = 2.*(gp[i] - ab[i]*cfact*od*ab2); | 
| 81 | 
+ | 
        } | 
| 82 | 
+ | 
        cnst = cfact*od*od*ab2 - DOT(gp,gp); | 
| 83 | 
+ | 
        /* | 
| 84 | 
+ | 
         * CONE:        sqcoef[0]*x*x + sqcoef[1]*y*y + sqcoef[2]*z*z | 
| 85 | 
+ | 
         *              + ctcoef[0]*x*y + ctcoef[1]*y*z + ctcoef[2]*z*x | 
| 86 | 
+ | 
         *              + licoef[0]*x + licoef[1]*y + licoef[2]*z + cnst == 0 | 
| 87 | 
+ | 
         */ | 
| 88 | 
+ | 
                                /* equation for conic section in plane */ | 
| 89 | 
+ | 
        gi[0] = hdwg0[gcl->gc.w]; | 
| 90 | 
+ | 
        gi[1] = hdwg1[gcl->gc.w]; | 
| 91 | 
+ | 
        wallpos = gcl->gc.w&1 ? gcl->hp->grid[gw] : 0; | 
| 92 | 
+ | 
        a = sqcoef[gi[0]];                                      /* x2 */ | 
| 93 | 
+ | 
        b = ctcoef[gi[0]];                                      /* xy */ | 
| 94 | 
+ | 
        c = sqcoef[gi[1]];                                      /* y2 */ | 
| 95 | 
+ | 
        d = ctcoef[gw]*wallpos + licoef[gi[0]];                 /* x */ | 
| 96 | 
+ | 
        e = ctcoef[gi[1]]*wallpos + licoef[gi[1]];              /* y */ | 
| 97 | 
+ | 
        f = wallpos*(wallpos*sqcoef[gw] + licoef[gw]) + cnst; | 
| 98 | 
+ | 
        for (i = 0; i < 2; i++) { | 
| 99 | 
+ | 
                if (i) {                /* swap x and y coefficients */ | 
| 100 | 
+ | 
                        register double t; | 
| 101 | 
+ | 
                        t = a; a = c; c = t; | 
| 102 | 
+ | 
                        t = d; d = e; e = t; | 
| 103 | 
+ | 
                } | 
| 104 | 
+ | 
                nex = 0;                /* check global extrema */ | 
| 105 | 
+ | 
                n = quadratic(root, a*(4.*a*c-b*b), 2.*a*(2.*c*d-b*e), | 
| 106 | 
+ | 
                                d*(c*d-b*e) + f*b*b); | 
| 107 | 
+ | 
                while (n-- > 0) { | 
| 108 | 
+ | 
                        if (gc->w>>1 == gi[i] && | 
| 109 | 
+ | 
                                        (gc->w&1) ^ root[n] < gp[gc->w>>1]) { | 
| 110 | 
+ | 
                                if (gc->w&1) | 
| 111 | 
+ | 
                                        gcl->gmin[i] = -FHUGE; | 
| 112 | 
+ | 
                                else | 
| 113 | 
+ | 
                                        gcl->gmax[i] = FHUGE; | 
| 114 | 
+ | 
                                nex++; | 
| 115 | 
+ | 
                                continue;               /* hyperbolic */ | 
| 116 | 
+ | 
                        } | 
| 117 | 
+ | 
                        if (tight) { | 
| 118 | 
+ | 
                                yex = (-2.*a*root[n] - d)/b; | 
| 119 | 
+ | 
                                if (yex < gcl->gc.i[1-i] || | 
| 120 | 
+ | 
                                                yex > gcl->gc.i[1-i]+1) | 
| 121 | 
+ | 
                                        continue;       /* outside cell */ | 
| 122 | 
+ | 
                        } | 
| 123 | 
+ | 
                        if (root[n] < gcl->gmin[i]) | 
| 124 | 
+ | 
                                gcl->gmin[i] = root[n]; | 
| 125 | 
+ | 
                        if (root[n] > gcl->gmax[i]) | 
| 126 | 
+ | 
                                gcl->gmax[i] = root[n]; | 
| 127 | 
+ | 
                        nex++; | 
| 128 | 
+ | 
                } | 
| 129 | 
+ | 
                                        /* check local extrema */ | 
| 130 | 
+ | 
                for (j = nex < 2 ? 2 : 0; j--; ) { | 
| 131 | 
+ | 
                        yex = gcl->gc.i[1-i] + j; | 
| 132 | 
+ | 
                        n = quadratic(root, a, b*yex+d, yex*(yex*c+e)+f); | 
| 133 | 
+ | 
                        while (n-- > 0) { | 
| 134 | 
+ | 
                                if (gc->w>>1 == gi[i] && | 
| 135 | 
+ | 
                                        (gc->w&1) ^ root[n] < gp[gc->w>>1]) | 
| 136 | 
+ | 
                                        continue; | 
| 137 | 
+ | 
                                if (root[n] < gcl->gmin[i]) | 
| 138 | 
+ | 
                                        gcl->gmin[i] = root[n]; | 
| 139 | 
+ | 
                                if (root[n] > gcl->gmax[i]) | 
| 140 | 
+ | 
                                        gcl->gmax[i] = root[n]; | 
| 141 | 
+ | 
                        } | 
| 142 | 
+ | 
                } | 
| 143 | 
+ | 
        } | 
| 144 | 
+ | 
} | 
| 145 | 
+ | 
 | 
| 146 | 
+ | 
 | 
| 147 | 
+ | 
static int | 
| 148 | 
+ | 
clipeyelim(rrng, gcl)           /* clip eye limits to grid cell */ | 
| 149 | 
+ | 
register short  rrng[2][2]; | 
| 150 | 
+ | 
register struct gclim   *gcl; | 
| 151 | 
+ | 
{ | 
| 152 | 
+ | 
        int     incell = 1; | 
| 153 | 
+ | 
        register int    i; | 
| 154 | 
+ | 
 | 
| 155 | 
+ | 
        for (i = 0; i < 2; i++) { | 
| 156 | 
+ | 
                if (gcl->gmin[i] < gcl->gc.i[i]) | 
| 157 | 
+ | 
                        gcl->gmin[i] = gcl->gc.i[i]; | 
| 158 | 
+ | 
                if (gcl->gmax[i] > gcl->gc.i[i]+1) | 
| 159 | 
+ | 
                        gcl->gmax[i] = gcl->gc.i[i]+1; | 
| 160 | 
+ | 
                if (gcl->gmax[i] > gcl->gmin[i]) { | 
| 161 | 
+ | 
                        rrng[i][0] = 256.*(gcl->gmin[i] - gcl->gc.i[i]) + | 
| 162 | 
+ | 
                                        (1.-FTINY); | 
| 163 | 
+ | 
                        rrng[i][1] = 256.*(gcl->gmax[i] - gcl->gc.i[i]) + | 
| 164 | 
+ | 
                                        (1.-FTINY) - rrng[i][0]; | 
| 165 | 
+ | 
                } else | 
| 166 | 
+ | 
                        rrng[i][0] = rrng[i][1] = 0; | 
| 167 | 
+ | 
                incell &= rrng[i][1] > 0; | 
| 168 | 
+ | 
        } | 
| 169 | 
+ | 
        return(incell); | 
| 170 | 
+ | 
} | 
| 171 | 
+ | 
 | 
| 172 | 
+ | 
 | 
| 173 | 
  | 
packrays(rod, p)                /* pack ray origins and directions */ | 
| 174 | 
  | 
register float  *rod; | 
| 175 | 
  | 
register PACKET *p; | 
| 176 | 
  | 
{ | 
| 177 | 
< | 
        short   packord[RPACKSIZ]; | 
| 178 | 
< | 
        float   packdc[RPACKSIZ]; | 
| 179 | 
< | 
        int     iterleft = 3*p->nr; | 
| 180 | 
< | 
        BYTE    rpos[2][2]; | 
| 181 | 
< | 
        FVECT   ro, rd, rp1; | 
| 177 | 
> | 
#if 0 | 
| 178 | 
> | 
        double  dist2sum = 0.; | 
| 179 | 
> | 
        FVECT   vt; | 
| 180 | 
> | 
#endif | 
| 181 | 
> | 
        int     nretries = p->nr + 2; | 
| 182 | 
> | 
        struct gclim    eyelim; | 
| 183 | 
> | 
        short   rrng0[2][2], rrng1[2][2]; | 
| 184 | 
> | 
        int     useyelim; | 
| 185 | 
  | 
        GCOORD  gc[2]; | 
| 186 | 
< | 
        double  d, dc, meandist; | 
| 187 | 
< | 
        int     i; | 
| 188 | 
< | 
        register int    ii; | 
| 186 | 
> | 
        FVECT   ro, rd; | 
| 187 | 
> | 
        double  d; | 
| 188 | 
> | 
        register int    i; | 
| 189 | 
  | 
 | 
| 190 | 
  | 
        if (!hdbcoord(gc, hdlist[p->hd], p->bi)) | 
| 191 | 
  | 
                error(CONSISTENCY, "bad beam index in packrays"); | 
| 192 | 
< | 
        for (i = 0, meandist = 0.; i < p->nr || meandist > myeye.rng+FTINY; ) { | 
| 193 | 
< | 
                rpos[0][0] = frandom() * 256.; | 
| 194 | 
< | 
                rpos[0][1] = frandom() * 256.; | 
| 195 | 
< | 
                rpos[1][0] = frandom() * 256.; | 
| 196 | 
< | 
                rpos[1][1] = frandom() * 256.; | 
| 197 | 
< | 
                d = hdray(ro, rd, hdlist[p->hd], gc, rpos); | 
| 198 | 
< | 
                if (myeye.rng > FTINY) {                /* check eyepoint */ | 
| 199 | 
< | 
                        register int    nexti; | 
| 200 | 
< | 
 | 
| 201 | 
< | 
                        VSUM(rp1, ro, rd, d); | 
| 202 | 
< | 
                        dc = sqrt(dist2line(myeye.vpt, ro, rp1)) / p->nr; | 
| 203 | 
< | 
                        if (i == p->nr) {               /* packet full */ | 
| 204 | 
< | 
                                nexti = packord[i-1]; | 
| 205 | 
< | 
                                if (!iterleft--) | 
| 206 | 
< | 
                                        break;          /* tried enough! */ | 
| 207 | 
< | 
                                if (dc >= packdc[nexti]) | 
| 208 | 
< | 
                                        continue;       /* worse than worst */ | 
| 209 | 
< | 
                                meandist -= packdc[nexti]; | 
| 210 | 
< | 
                        } else | 
| 211 | 
< | 
                                nexti = i++; | 
| 212 | 
< | 
                        meandist += packdc[nexti] = dc; /* new distance */ | 
| 213 | 
< | 
                        for (ii = i; --ii; ) {          /* insertion sort */ | 
| 214 | 
< | 
                                if (dc > packdc[packord[ii-1]]) | 
| 215 | 
< | 
                                        break; | 
| 216 | 
< | 
                                packord[ii] = packord[ii-1]; | 
| 192 | 
> | 
        if ((useyelim = myeye.rng > FTINY)) { | 
| 193 | 
> | 
                initeyelim(&eyelim, hdlist[p->hd], gc); | 
| 194 | 
> | 
                groweyelim(&eyelim, gc+1, 0., 0., 0); | 
| 195 | 
> | 
                groweyelim(&eyelim, gc+1, 1., 1., 0); | 
| 196 | 
> | 
                useyelim = clipeyelim(rrng0, &eyelim); | 
| 197 | 
> | 
#ifdef DEBUG | 
| 198 | 
> | 
                if (!useyelim) | 
| 199 | 
> | 
                        error(WARNING, "no eye overlap in packrays"); | 
| 200 | 
> | 
#endif | 
| 201 | 
> | 
        } | 
| 202 | 
> | 
        for (i = 0; i < p->nr; i++) { | 
| 203 | 
> | 
        retry: | 
| 204 | 
> | 
                if (useyelim) { | 
| 205 | 
> | 
                        initeyelim(&eyelim, NULL, gc+1); | 
| 206 | 
> | 
                        p->ra[i].r[0][0] = (int)(frandom()*rrng0[0][1]) | 
| 207 | 
> | 
                                                + rrng0[0][0]; | 
| 208 | 
> | 
                        p->ra[i].r[0][1] = (int)(frandom()*rrng0[1][1]) | 
| 209 | 
> | 
                                                + rrng0[1][0]; | 
| 210 | 
> | 
                        groweyelim(&eyelim, gc, | 
| 211 | 
> | 
                                        (1./256.)*(p->ra[i].r[0][0]+.5), | 
| 212 | 
> | 
                                        (1./256.)*(p->ra[i].r[0][1]+.5), 1); | 
| 213 | 
> | 
                        if (!clipeyelim(rrng1, &eyelim)) { | 
| 214 | 
> | 
                                useyelim = nretries-- > 0; | 
| 215 | 
> | 
#ifdef DEBUG | 
| 216 | 
> | 
                                if (!useyelim) | 
| 217 | 
> | 
                                        error(WARNING, | 
| 218 | 
> | 
                                        "exceeded retry limit in packrays"); | 
| 219 | 
> | 
#endif | 
| 220 | 
> | 
                                goto retry; | 
| 221 | 
  | 
                        } | 
| 222 | 
< | 
                        packord[ii] = nexti; | 
| 223 | 
< | 
                        ii = nexti;                     /* put it here */ | 
| 224 | 
< | 
                } else | 
| 225 | 
< | 
                        ii = i++; | 
| 222 | 
> | 
                        p->ra[i].r[1][0] = (int)(frandom()*rrng1[0][1]) | 
| 223 | 
> | 
                                                + rrng1[0][0]; | 
| 224 | 
> | 
                        p->ra[i].r[1][1] = (int)(frandom()*rrng1[1][1]) | 
| 225 | 
> | 
                                                + rrng1[1][0]; | 
| 226 | 
> | 
                } else { | 
| 227 | 
> | 
                        p->ra[i].r[0][0] = frandom() * 256.; | 
| 228 | 
> | 
                        p->ra[i].r[0][1] = frandom() * 256.; | 
| 229 | 
> | 
                        p->ra[i].r[1][0] = frandom() * 256.; | 
| 230 | 
> | 
                        p->ra[i].r[1][1] = frandom() * 256.; | 
| 231 | 
> | 
                } | 
| 232 | 
> | 
                d = hdray(ro, rd, hdlist[p->hd], gc, p->ra[i].r); | 
| 233 | 
> | 
#if 0 | 
| 234 | 
> | 
                VSUM(vt, ro, rd, d); | 
| 235 | 
> | 
                dist2sum += dist2line(myeye.vpt, ro, vt); | 
| 236 | 
> | 
#endif | 
| 237 | 
  | 
                if (p->offset != NULL) { | 
| 238 | 
  | 
                        if (!vdef(OBSTRUCTIONS)) | 
| 239 | 
  | 
                                d *= frandom();         /* random offset */ | 
| 240 | 
  | 
                        VSUM(ro, ro, rd, d);            /* advance ray */ | 
| 241 | 
< | 
                        p->offset[ii] = d; | 
| 241 | 
> | 
                        p->offset[i] = d; | 
| 242 | 
  | 
                } | 
| 243 | 
< | 
                p->ra[ii].r[0][0] = rpos[0][0]; | 
| 244 | 
< | 
                p->ra[ii].r[0][1] = rpos[0][1]; | 
| 245 | 
< | 
                p->ra[ii].r[1][0] = rpos[1][0]; | 
| 246 | 
< | 
                p->ra[ii].r[1][1] = rpos[1][1]; | 
| 75 | 
< | 
                VCOPY(rod+6*ii, ro); | 
| 76 | 
< | 
                VCOPY(rod+6*ii+3, rd); | 
| 243 | 
> | 
                VCOPY(rod, ro); | 
| 244 | 
> | 
                rod += 3; | 
| 245 | 
> | 
                VCOPY(rod, rd); | 
| 246 | 
> | 
                rod += 3; | 
| 247 | 
  | 
        } | 
| 248 | 
< | 
#ifdef DEBUG | 
| 249 | 
< | 
        fprintf(stderr, "%f mean distance for target %f (%d iterations)\n", | 
| 250 | 
< | 
                        meandist, myeye.rng, 3*p->nr - iterleft); | 
| 248 | 
> | 
#if 0 | 
| 249 | 
> | 
        fprintf(stderr, "%f RMS (%d retries)\t", sqrt(dist2sum/p->nr), | 
| 250 | 
> | 
                        p->nr + 2 - nretries); | 
| 251 | 
  | 
#endif | 
| 252 | 
  | 
} | 
| 253 | 
  | 
 | 
| 277 | 
  | 
        int     status; | 
| 278 | 
  | 
                                        /* already closed? */ | 
| 279 | 
  | 
        if (!nprocs) | 
| 280 | 
< | 
                return; | 
| 280 | 
> | 
                return(0); | 
| 281 | 
  | 
                                        /* flush beam queue */ | 
| 282 | 
  | 
        done_packets(flush_queue()); | 
| 283 | 
  | 
                                        /* sync holodeck */ | 
| 348 | 
  | 
                        pippt = NULL; | 
| 349 | 
  | 
                } | 
| 350 | 
  | 
        if (pippt != NULL) | 
| 351 | 
< | 
                strcpy(pippt, "> /dev/null");   /* nothing to match */ | 
| 351 | 
> | 
                strcpy(pippt, "> " NULL_DEVICE);        /* nothing to match */ | 
| 352 | 
  | 
        else | 
| 353 | 
  | 
                sprintf(cp, ")[ \t]*=' > %s", tf2); | 
| 354 | 
  | 
#ifdef DEBUG |