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

Comparing ray/src/hd/rholo2.c (file contents):
Revision 3.4 by gregl, Mon Dec 1 16:34:36 1997 UTC vs.
Revision 3.19 by gwlarson, Tue Dec 1 15:47:05 1998 UTC

# Line 1 | Line 1
1 < /* Copyright (c) 1997 Silicon Graphics, Inc. */
1 > /* Copyright (c) 1998 Silicon Graphics, Inc. */
2  
3   #ifndef lint
4   static char SCCSid[] = "$SunId$ SGI";
# Line 9 | Line 9 | static char SCCSid[] = "$SunId$ SGI";
9   */
10  
11   #include "rholo.h"
12 + #include "paths.h"
13   #include "random.h"
14  
15  
16 + VIEWPOINT       myeye;          /* target view position */
17 +
18 + struct gclim {
19 +        HOLO    *hp;                    /* holodeck pointer */
20 +        GCOORD  gc;                     /* grid cell */
21 +        FVECT   egp;                    /* eye grid point */
22 +        double  erg[2];                 /* eye range in wall grid coords */
23 +        double  gmin[2], gmax[2];       /* grid coordinate limits */
24 + };                              /* a grid coordinate range */
25 +
26 +
27 + static
28 + initeyelim(gcl, hd, gc)         /* initialize grid coordinate limits */
29 + register struct gclim   *gcl;
30 + int     hd;
31 + GCOORD  *gc;
32 + {
33 +        register FLOAT  *v;
34 +        register int    i;
35 +
36 +        gcl->hp = hdlist[hd];
37 +        copystruct(&gcl->gc, gc);
38 +        hdgrid(gcl->egp, gcl->hp, myeye.vpt);
39 +        for (i = 0; i < 2; i++) {
40 +                v = gcl->hp->wg[((gcl->gc.w>>1)+i+1)%3];
41 +                gcl->erg[i] = myeye.rng * VLEN(v);
42 +                gcl->gmin[i] = FHUGE; gcl->gmax[i] = -FHUGE;
43 +        }
44 + }
45 +
46 +
47 + static
48 + groweyelim(gcl, gp)             /* grow grid limits about eye point */
49 + register struct gclim   *gcl;
50 + FVECT   gp;
51 + {
52 +        FVECT   ab;
53 +        double  l2, d, mult, wg;
54 +        register int    i, g;
55 +
56 +        VSUB(ab, gcl->egp, gp);
57 +        l2 = DOT(ab,ab);
58 +        if (l2 <= gcl->erg[0]*gcl->erg[1]) {
59 +                gcl->gmin[0] = gcl->gmin[1] = -FHUGE;
60 +                gcl->gmax[0] = gcl->gmax[1] = FHUGE;
61 +                return;
62 +        }
63 +        mult = gp[g = gcl->gc.w>>1];
64 +        if (gcl->gc.w&1)
65 +                mult -= gcl->hp->grid[g];
66 +        if (ab[g]*ab[g] > gcl->erg[0]*gcl->erg[1])
67 +                mult /= -ab[g];
68 +        else if (fabs(ab[hdwg0[gcl->gc.w]]) > fabs(ab[hdwg1[gcl->gc.w]]))
69 +                mult = (gcl->gc.i[0] + .5 - gp[hdwg0[gcl->gc.w]]) /
70 +                                ab[hdwg0[gcl->gc.w]];
71 +        else
72 +                mult = (gcl->gc.i[1] + .5 - gp[hdwg1[gcl->gc.w]]) /
73 +                                ab[hdwg1[gcl->gc.w]];
74 +        for (i = 0; i < 2; i++) {
75 +                g = ((gcl->gc.w>>1)+i+1)%3;
76 +                wg = gp[g] + mult*ab[g];
77 +                d = mult*gcl->erg[i];
78 +                if (d < 0.) d = -d;
79 +                if (wg - d < gcl->gmin[i])
80 +                        gcl->gmin[i] = wg - d;
81 +                if (wg + d > gcl->gmax[i])
82 +                        gcl->gmax[i] = wg + d;
83 +        }
84 + }
85 +
86 +
87 + static int
88 + clipeyelim(rrng, gcl)           /* clip eye limits to grid cell */
89 + register short  rrng[2][2];
90 + register struct gclim   *gcl;
91 + {
92 +        int     incell = 1;
93 +        register int    i;
94 +
95 +        for (i = 0; i < 2; i++) {
96 +                if (gcl->gmin[i] < gcl->gc.i[i])
97 +                        gcl->gmin[i] = gcl->gc.i[i];
98 +                if (gcl->gmax[i] > gcl->gc.i[i]+1)
99 +                        gcl->gmax[i] = gcl->gc.i[i]+1;
100 +                if ((incell &= gcl->gmax[i] > gcl->gmin[i])) {
101 +                        rrng[i][0] = 256.*(gcl->gmin[i] - gcl->gc.i[i]) +
102 +                                        (1.-FTINY);
103 +                        rrng[i][1] = 256.*(gcl->gmax[i] - gcl->gc.i[i]) +
104 +                                        (1.-FTINY) - rrng[i][0];
105 +                        incell &= rrng[i][1] > 0;
106 +                }
107 +        }
108 +        return(incell);
109 + }
110 +
111 +
112   packrays(rod, p)                /* pack ray origins and directions */
113   register float  *rod;
114   register PACKET *p;
115   {
116 <        static FVECT    ro, rd;
116 > #define gp      ro
117 > #ifdef DEBUG
118 >        double dist2sum = 0.;
119 >        FVECT   vt;
120 > #endif
121 >        int     nretries = p->nr + 2;
122 >        struct gclim    eyelim;
123 >        short   rrng0[2][2], rrng1[2][2];
124 >        int     useyelim;
125          GCOORD  gc[2];
126 <        int     ila[2], hsh;
127 <        double  d, sl[4];
126 >        FVECT   ro, rd;
127 >        double  d;
128          register int    i;
129  
130          if (!hdbcoord(gc, hdlist[p->hd], p->bi))
131                  error(CONSISTENCY, "bad beam index in packrays");
132 <        ila[0] = p->hd; ila[1] = p->bi;
133 <        hsh = ilhash(ila,2) + p->nc;
132 >        if ((useyelim = myeye.rng > FTINY)) {
133 >                initeyelim(&eyelim, p->hd, gc);
134 >                gp[gc[1].w>>1] = gc[1].w&1 ?
135 >                                hdlist[p->hd]->grid[gc[1].w>>1] : 0;
136 >                gp[hdwg0[gc[1].w]] = gc[1].i[0];
137 >                gp[hdwg1[gc[1].w]] = gc[1].i[1];
138 >                groweyelim(&eyelim, gp);
139 >                gp[hdwg0[gc[1].w]]++;
140 >                gp[hdwg1[gc[1].w]]++;
141 >                groweyelim(&eyelim, gp);
142 >                useyelim &= clipeyelim(rrng0, &eyelim);
143 >        }
144          for (i = 0; i < p->nr; i++) {
145 <                multisamp(sl, 4, urand(hsh+i));
146 <                p->ra[i].r[0][0] = sl[0] * 256.;
147 <                p->ra[i].r[0][1] = sl[1] * 256.;
148 <                p->ra[i].r[1][0] = sl[2] * 256.;
149 <                p->ra[i].r[1][1] = sl[3] * 256.;
145 >        retry:
146 >                if (useyelim) {
147 >                        p->ra[i].r[0][0] = (int)(frandom()*rrng0[0][1])
148 >                                                + rrng0[0][0];
149 >                        p->ra[i].r[0][1] = (int)(frandom()*rrng0[1][1])
150 >                                                + rrng0[1][0];
151 >                        initeyelim(&eyelim, p->hd, gc+1);
152 >                        gp[gc[0].w>>1] = gc[0].w&1 ?
153 >                                        hdlist[p->hd]->grid[gc[0].w>>1] : 0;
154 >                        gp[hdwg0[gc[0].w]] = gc[0].i[0] +
155 >                                        (1./256.)*(p->ra[i].r[0][0]+.5);
156 >                        gp[hdwg1[gc[0].w]] = gc[0].i[1] +
157 >                                        (1./256.)*(p->ra[i].r[0][1]+.5);
158 >                        groweyelim(&eyelim, gp);
159 >                        if (!clipeyelim(rrng1, &eyelim)) {
160 >                                useyelim &= nretries-- > 0;
161 > #ifdef DEBUG
162 >                                if (!useyelim)
163 >                                        error(WARNING, "exceeded retry limit in packrays");
164 > #endif
165 >                                goto retry;
166 >                        }
167 >                        p->ra[i].r[1][0] = (int)(frandom()*rrng1[0][1])
168 >                                                + rrng1[0][0];
169 >                        p->ra[i].r[1][1] = (int)(frandom()*rrng1[1][1])
170 >                                                + rrng1[1][0];
171 >                } else {
172 >                        p->ra[i].r[0][0] = frandom() * 256.;
173 >                        p->ra[i].r[0][1] = frandom() * 256.;
174 >                        p->ra[i].r[1][0] = frandom() * 256.;
175 >                        p->ra[i].r[1][1] = frandom() * 256.;
176 >                }
177                  d = hdray(ro, rd, hdlist[p->hd], gc, p->ra[i].r);
178 + #ifdef DEBUG
179 +                VSUM(vt, ro, rd, d);
180 +                dist2sum += dist2line(myeye.vpt, ro, vt);
181 + #endif
182                  if (p->offset != NULL) {
183 <                        VSUM(ro, ro, rd, d);            /* exterior only */
183 >                        if (!vdef(OBSTRUCTIONS))
184 >                                d *= frandom();         /* random offset */
185 >                        VSUM(ro, ro, rd, d);            /* advance ray */
186                          p->offset[i] = d;
187                  }
188                  VCOPY(rod, ro);
# Line 42 | Line 190 | register PACKET        *p;
190                  VCOPY(rod, rd);
191                  rod += 3;
192          }
193 + #ifdef DEBUG
194 +        fprintf(stderr, "RMS distance = %f\n", sqrt(dist2sum/p->nr));
195 + #endif
196 + #undef gp
197   }
198  
199  
# Line 61 | Line 213 | register float *rvl;
213                  rvl += 4;
214          }
215          p->nc += p->nr;
216 + }
217 +
218 +
219 + int
220 + done_rtrace()                   /* clean up and close rtrace calculation */
221 + {
222 +        int     status;
223 +                                        /* already closed? */
224 +        if (!nprocs)
225 +                return;
226 +                                        /* flush beam queue */
227 +        done_packets(flush_queue());
228 +                                        /* sync holodeck */
229 +        hdsync(NULL, 1);
230 +                                        /* close rtrace */
231 +        if ((status = end_rtrace()))
232 +                error(WARNING, "bad exit status from rtrace");
233 +        if (vdef(REPORT)) {             /* report time */
234 +                eputs("rtrace process closed\n");
235 +                report(0);
236 +        }
237 +        return(status);                 /* return status */
238 + }
239 +
240 +
241 + new_rtrace()                    /* restart rtrace calculation */
242 + {
243 +        char    combuf[128];
244 +
245 +        if (nprocs > 0)                 /* already running? */
246 +                return;
247 +        starttime = time(NULL);         /* reset start time and counts */
248 +        npacksdone = nraysdone = 0L;
249 +        if (vdef(TIME))                 /* reset end time */
250 +                endtime = starttime + vflt(TIME)*3600. + .5;
251 +        if (vdef(RIF)) {                /* rerun rad to update octree */
252 +                sprintf(combuf, "rad -v 0 -s -w %s", vval(RIF));
253 +                if (system(combuf))
254 +                        error(WARNING, "error running rad");
255 +        }
256 +        if (start_rtrace() < 1)         /* start rtrace */
257 +                error(WARNING, "cannot restart rtrace");
258 +        else if (vdef(REPORT)) {
259 +                eputs("rtrace process restarted\n");
260 +                report(0);
261 +        }
262 + }
263 +
264 +
265 + getradfile()                    /* run rad and get needed variables */
266 + {
267 +        static short    mvar[] = {OCTREE,EYESEP,-1};
268 +        static char     tf1[] = TEMPLATE;
269 +        char    tf2[64];
270 +        char    combuf[256];
271 +        char    *pippt;
272 +        register int    i;
273 +        register char   *cp;
274 +                                        /* check if rad file specified */
275 +        if (!vdef(RIF))
276 +                return(0);
277 +                                        /* create rad command */
278 +        mktemp(tf1);
279 +        sprintf(tf2, "%s.rif", tf1);
280 +        sprintf(combuf,
281 +                "rad -v 0 -s -e -w %s OPTFILE=%s | egrep '^[ \t]*(NOMATCH",
282 +                        vval(RIF), tf1);
283 +        cp = combuf;
284 +        while (*cp){
285 +                if (*cp == '|') pippt = cp;
286 +                cp++;
287 +        }                               /* match unset variables */
288 +        for (i = 0; mvar[i] >= 0; i++)
289 +                if (!vdef(mvar[i])) {
290 +                        *cp++ = '|';
291 +                        strcpy(cp, vnam(mvar[i]));
292 +                        while (*cp) cp++;
293 +                        pippt = NULL;
294 +                }
295 +        if (pippt != NULL)
296 +                strcpy(pippt, "> /dev/null");   /* nothing to match */
297 +        else
298 +                sprintf(cp, ")[ \t]*=' > %s", tf2);
299 + #ifdef DEBUG
300 +        wputs(combuf); wputs("\n");
301 + #endif
302 +        system(combuf);                         /* ignore exit code */
303 +        if (pippt == NULL) {
304 +                loadvars(tf2);                  /* load variables */
305 +                unlink(tf2);
306 +        }
307 +        rtargc += wordfile(rtargv+rtargc, tf1); /* get rtrace options */
308 +        unlink(tf1);                    /* clean up */
309 +        return(1);
310 + }
311 +
312 +
313 + report(t)                       /* report progress so far */
314 + time_t  t;
315 + {
316 +        static time_t   seconds2go = 1000000;
317 +
318 +        if (t == 0L)
319 +                t = time(NULL);
320 +        sprintf(errmsg, "%ld packets (%ld rays) done after %.2f hours\n",
321 +                        npacksdone, nraysdone, (t-starttime)/3600.);
322 +        eputs(errmsg);
323 +        if (seconds2go == 1000000)
324 +                seconds2go = vdef(REPORT) ? (long)(vflt(REPORT)*60. + .5) : 0L;
325 +        if (seconds2go)
326 +                reporttime = t + seconds2go;
327   }

Diff Legend

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