ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/hd/holo.c
Revision: 3.2
Committed: Mon Nov 3 11:03:23 1997 UTC (26 years, 5 months ago) by gregl
Content type: text/plain
Branch: MAIN
Changes since 3.1: +30 -43 lines
Log Message:
changed to indexed array to avoid expensive modulo evaluations

File Contents

# User Rev Content
1 gregl 3.1 /* Copyright (c) 1997 Silicon Graphics, Inc. */
2    
3     #ifndef lint
4     static char SCCSid[] = "$SunId$ SGI";
5     #endif
6    
7     /*
8     * Routines for converting holodeck coordinates, etc.
9     *
10     * 10/22/97 GWLarson
11     */
12    
13     #include "holo.h"
14    
15     float hd_depthmap[DCINF-DCLIN];
16    
17     static double logstep;
18    
19 gregl 3.2 static int wg0[6] = {1,1,2,2,0,0};
20     static int wg1[6] = {2,2,0,0,1,1};
21 gregl 3.1
22 gregl 3.2
23 gregl 3.1 hdcompgrid(hp) /* compute derived grid vector and index */
24     register HOLO *hp;
25     {
26     FVECT AxB;
27     double d;
28     register FLOAT *v;
29     register int i, j;
30     /* initialize depth map */
31     if (hd_depthmap[0] < 1.) {
32     d = 1. + .5/DCLIN;
33     for (i = 0; i < DCINF-DCLIN; i++) {
34     hd_depthmap[i] = d;
35     d *= 1. + 1./DCLIN;
36     }
37     logstep = log(1. + 1./DCLIN);
38     }
39     /* compute grid coordinate vectors */
40     for (i = 0; i < 3; i++) {
41     fcross(AxB, hp->xv[(i+1)%3], v=hp->xv[(i+2)%3]);
42     VCOPY(hp->wn[i], AxB);
43     if (normalize(hp->wn[i]) == 0.)
44     error(USER, "degenerate holodeck section");
45     hp->wo[i<<1] = DOT(hp->wn[i],hp->orig);
46     hp->wo[i<<1|1] = hp->wo[i<<1] + DOT(hp->wn[i],hp->xv[i]);
47     fcross(hp->gv[i][0], v, AxB);
48     d = DOT(v,v) / DOT(hp->gv[i][0],hp->gv[i][0]) *
49     hp->grid[(i+1)%3];
50     for (j = 0; j < 3; j++)
51     hp->gv[i][0][j] *= d;
52     fcross(hp->gv[i][1], AxB, v=hp->xv[(i+1)%3]);
53     d = DOT(v,v) / DOT(hp->gv[i][1],hp->gv[i][1]) *
54     hp->grid[(i+2)%3];
55     for (j = 0; j < 3; j++)
56     hp->gv[i][1][j] *= d;
57     }
58     /* compute linear depth range */
59     hp->tlin = VLEN(hp->xv[0]) + VLEN(hp->xv[1]) + VLEN(hp->xv[2]);
60     /* compute wall super-indices from grid */
61     hp->wi[0] = 1; /**** index values begin at 1 ****/
62     for (i = 1; i < 6; i++) {
63     hp->wi[i] = 0;
64     for (j = i; j < 6; j++)
65 gregl 3.2 hp->wi[i] += hp->grid[wg0[j]] * hp->grid[wg1[j]];
66     hp->wi[i] *= hp->grid[wg0[i-1]] * hp->grid[wg1[i-1]];
67 gregl 3.1 hp->wi[i] += hp->wi[i-1];
68     }
69     }
70    
71    
72     HOLO *
73     hdalloc(hproto) /* allocate and set holodeck section based on grid */
74     HDGRID *hproto;
75     {
76     HOLO hdhead;
77     register HOLO *hp;
78     int n;
79     /* copy grid to temporary header */
80     bcopy((char *)hproto, (char *)&hdhead, sizeof(HDGRID));
81     /* compute grid vectors and sizes */
82     hdcompgrid(&hdhead);
83     /* allocate header with directory */
84     n = sizeof(HOLO)+nbeams(&hdhead)*sizeof(BEAMI);
85     if ((hp = (HOLO *)malloc(n)) == NULL)
86     return(NULL);
87     /* copy header information */
88     copystruct(hp, &hdhead);
89     /* allocate and clear beam list */
90     hp->bl = (BEAM **)malloc((nbeams(hp)+1)*sizeof(BEAM *)+sizeof(BEAM));
91     if (hp->bl == NULL) {
92     free((char *)hp);
93     return(NULL);
94     }
95     bzero((char *)hp->bl, (nbeams(hp)+1)*sizeof(BEAM *)+sizeof(BEAM));
96     hp->bl[0] = (BEAM *)(hp->bl+nbeams(hp)+1); /* set blglob(hp) */
97     hp->fd = -1;
98     hp->dirty = 0;
99     hp->priv = NULL;
100     /* clear beam directory */
101     bzero((char *)hp->bi, (nbeams(hp)+1)*sizeof(BEAMI));
102     return(hp); /* all is well */
103     }
104    
105    
106     hdbcoord(gc, hp, i) /* compute beam coordinates from index */
107     BCOORD gc; /* returned */
108     register HOLO *hp;
109     register int i;
110     {
111     register int j, n;
112     int n2, reverse;
113     BCOORD g2;
114     /* check range */
115     if (i < 1 | i > nbeams(hp))
116     return(0);
117     if (reverse = i >= hp->wi[5])
118     i -= hp->wi[5] - 1;
119     for (j = 0; j < 5; j++) /* find w0 */
120     if (hp->wi[j+1] > i)
121     break;
122     i -= hp->wi[gc[0].w=j];
123     /* find w1 */
124 gregl 3.2 n2 = hp->grid[wg0[j]] * hp->grid[wg1[j]];
125 gregl 3.1 while (++j < 5) {
126 gregl 3.2 n = n2 * hp->grid[wg0[j]] * hp->grid[wg1[j]];
127 gregl 3.1 if (n > i)
128     break;
129     i -= n;
130     }
131     gc[1].w = j;
132     /* find position on w0 */
133 gregl 3.2 n2 = hp->grid[wg0[j]] * hp->grid[wg1[j]];
134 gregl 3.1 n = i / n2;
135 gregl 3.2 gc[0].i[1] = n / hp->grid[wg0[gc[0].w]];
136     gc[0].i[0] = n - gc[0].i[1]*hp->grid[wg0[gc[0].w]];
137 gregl 3.1 i -= n*n2;
138     /* find position on w1 */
139 gregl 3.2 gc[1].i[1] = i / hp->grid[wg0[gc[1].w]];
140     gc[1].i[0] = i - gc[1].i[1]*hp->grid[wg0[gc[1].w]];
141 gregl 3.1 if (reverse) {
142     copystruct(g2, gc+1);
143     copystruct(gc+1, gc);
144     copystruct(gc, g2);
145     }
146     return(1); /* we're done */
147     }
148    
149    
150     int
151     hdbindex(hp, gc) /* compute index from beam coordinates */
152     register HOLO *hp;
153     register BCOORD gc;
154     {
155     BCOORD g2;
156     int reverse;
157     register int i, j;
158     /* check ordering and limits */
159     if (reverse = gc[0].w > gc[1].w) {
160     copystruct(g2, gc+1);
161     copystruct(g2+1, gc);
162     gc = g2;
163     } else if (gc[0].w == gc[1].w)
164     return(0);
165     if (gc[0].w < 0 | gc[1].w > 5)
166     return(0);
167     i = 0; /* compute index */
168     for (j = gc[0].w+1; j < gc[1].w; j++)
169 gregl 3.2 i += hp->grid[wg0[j]] * hp->grid[wg1[j]];
170     i *= hp->grid[wg0[gc[0].w]] * hp->grid[wg1[gc[0].w]];
171 gregl 3.1 i += hp->wi[gc[0].w];
172 gregl 3.2 i += (hp->grid[wg0[gc[0].w]]*gc[0].i[1] + gc[0].i[0]) *
173     hp->grid[wg0[gc[1].w]] * hp->grid[wg1[gc[1].w]] ;
174     i += hp->grid[wg0[gc[1].w]]*gc[1].i[1] + gc[1].i[0];
175 gregl 3.1 if (reverse)
176     i += hp->wi[5] - 1;
177     return(i);
178     }
179    
180    
181     hdlseg(lseg, hp, i) /* compute line segment for beam */
182 gregl 3.2 register int lseg[2][3];
183 gregl 3.1 register HOLO *hp;
184     int i;
185     {
186     BCOORD gc;
187     register int k;
188    
189 gregl 3.2 if (!hdbcoord(gc, hp, i)) /* compute grid coordinates */
190 gregl 3.1 return(0);
191 gregl 3.2 for (k = 0; k < 2; k++) { /* compute end points */
192     lseg[k][gc[k].w>>1] = gc[k].w&1 ? hp->grid[gc[k].w>>1]-1 : 0 ;
193     lseg[k][wg0[gc[k].w]] = gc[k].i[0];
194     lseg[k][wg1[gc[k].w]] = gc[k].i[1];
195     }
196 gregl 3.1 return(1);
197     }
198    
199    
200     unsigned
201     hdcode(hp, d) /* compute depth code for d */
202     HOLO *hp;
203     double d;
204     {
205     double tl = hp->tlin;
206     register unsigned c;
207    
208     if (d <= 0.)
209     return(0);
210     if (d >= .99*FHUGE)
211     return(DCINF);
212     if (d < tl)
213     return((unsigned)(d*DCLIN/tl));
214     c = (unsigned)(log(d/tl)/logstep) + DCLIN;
215     return(c > DCINF ? DCINF : c);
216     }
217    
218    
219     double
220     hdray(ro, rd, hp, gc, r) /* compute ray within a beam */
221     FVECT ro, rd; /* returned */
222     register HOLO *hp;
223     register BCOORD gc;
224     BYTE r[2][2];
225     {
226     FVECT p[2];
227     register int i;
228     register FLOAT *v;
229     double d;
230     /* compute entry and exit points */
231     for (i = 0; i < 2; i++) {
232     VCOPY(p[i], hp->orig);
233     if (gc[i].w & 1) {
234     v = hp->xv[gc[i].w>>1];
235     p[i][0] += *v++; p[i][1] += *v++; p[i][2] += *v;
236     }
237     d = ( gc[i].i[0] + (1./256.)*(r[i][0]+.5) ) /
238 gregl 3.2 hp->grid[wg0[gc[i].w]];
239     v = hp->xv[wg0[gc[i].w]];
240 gregl 3.1 p[i][0] += d * *v++; p[i][1] += d * *v++; p[i][2] += d * *v;
241     d = (gc[i].i[1] + (1./256.)*(r[i][1]+.5)) /
242 gregl 3.2 hp->grid[wg1[gc[i].w]];
243     v = hp->xv[wg1[gc[i].w]];
244 gregl 3.1 p[i][0] += d * *v++; p[i][1] += d * *v++; p[i][2] += d * *v;
245     }
246     VCOPY(ro, p[0]); /* assign ray origin and direction */
247     rd[0] = p[1][0] - p[0][0];
248     rd[1] = p[1][1] - p[0][1];
249     rd[2] = p[1][2] - p[0][2];
250     return(normalize(rd)); /* return maximum inside distance */
251     }
252    
253    
254     double
255     hdinter(gc, r, hp, ro, rd) /* compute ray intersection with section */
256     register BCOORD gc; /* returned */
257     BYTE r[2][2]; /* returned */
258     register HOLO *hp;
259     FVECT ro, rd; /* rd should be normalized */
260     {
261     FVECT p[2], vt;
262     double d, t0, t1, d0, d1;
263     register FLOAT *v;
264     register int i;
265     /* first, intersect walls */
266     gc[0].w = gc[1].w = -1;
267     t0 = -FHUGE; t1 = FHUGE;
268     for (i = 0; i < 3; i++) { /* for each wall pair */
269     d = -DOT(rd, hp->wn[i]); /* plane distance */
270     if (d <= FTINY && d >= -FTINY) /* check for parallel */
271     continue;
272     d1 = DOT(ro, hp->wn[i]); /* ray distances */
273     d0 = (d1 - hp->wo[i<<1]) / d;
274     d1 = (d1 - hp->wo[i<<1|1]) / d;
275     if (d0 < d1) { /* check against best */
276     if (d0 > t0) {
277     t0 = d0;
278     gc[0].w = i<<1;
279     }
280     if (d1 < t1) {
281     t1 = d1;
282     gc[1].w = i<<1 | 1;
283     }
284     } else {
285     if (d1 > t0) {
286     t0 = d1;
287     gc[0].w = i<<1 | 1;
288     }
289     if (d0 < t1) {
290     t1 = d0;
291     gc[1].w = i<<1;
292     }
293     }
294     }
295     if (gc[0].w < 0 | gc[1].w < 0) /* paranoid check */
296     return(FHUGE);
297     /* compute intersections */
298     for (i = 0; i < 3; i++) {
299     p[0][i] = ro[i] + rd[i]*t0;
300     p[1][i] = ro[i] + rd[i]*t1;
301     }
302     /* now, compute grid coordinates */
303     for (i = 0; i < 2; i++) {
304     vt[0] = p[i][0] - hp->orig[0];
305     vt[1] = p[i][1] - hp->orig[1];
306     vt[2] = p[i][2] - hp->orig[2];
307     if (gc[i].w & 1) {
308     v = hp->xv[gc[i].w>>1];
309     vt[0] -= *v++; vt[1] -= *v++; vt[2] -= *v;
310     }
311     v = hp->gv[gc[i].w>>1][0];
312     d = DOT(vt, v);
313 gregl 3.2 if (d < 0. || (gc[i].i[0] = d) >= hp->grid[wg0[gc[i].w]])
314 gregl 3.1 return(FHUGE); /* outside wall */
315     r[i][0] = 256. * (d - gc[i].i[0]);
316     v = hp->gv[gc[i].w>>1][1];
317     d = DOT(vt, v);
318 gregl 3.2 if (d < 0. || (gc[i].i[1] = d) >= hp->grid[wg1[gc[i].w]])
319 gregl 3.1 return(FHUGE); /* outside wall */
320     r[i][1] = 256. * (d - gc[i].i[1]);
321     }
322     /* return distance from entry point */
323     vt[0] = ro[0] - p[0][0];
324     vt[1] = ro[1] - p[0][1];
325     vt[2] = ro[2] - p[0][2];
326     return(DOT(vt,rd));
327     }