ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/hd/rhdisp3.c
Revision: 3.16
Committed: Tue Jun 8 19:48:30 2004 UTC (19 years, 9 months ago) by greg
Content type: text/plain
Branch: MAIN
CVS Tags: rad4R2P2, rad5R0, rad5R1, rad3R7P2, rad3R7P1, rad4R2, rad4R1, rad4R0, rad3R6, rad3R6P1, rad3R8, rad3R9, rad4R2P1
Changes since 3.15: +1 -2 lines
Log Message:
Removed redundant #include's and fixed ordering on some headers

File Contents

# User Rev Content
1 gregl 3.1 #ifndef lint
2 greg 3.16 static const char RCSid[] = "$Id: rhdisp3.c,v 3.15 2004/01/01 11:21:55 schorsch Exp $";
3 gregl 3.1 #endif
4     /*
5 gregl 3.2 * Holodeck beam support for display process
6 gregl 3.1 */
7    
8     #include "rholo.h"
9     #include "rhdisp.h"
10    
11 gregl 3.4 struct cellist {
12     GCOORD *cl;
13     int n;
14     };
15 gregl 3.1
16 gregl 3.4
17 gregl 3.1 int
18 gwlarson 3.10 npixels(vp, hr, vr, hp, bi) /* compute appropriate nrays to evaluate */
19 gregl 3.5 register VIEW *vp;
20 gregl 3.1 int hr, vr;
21     HOLO *hp;
22     int bi;
23     {
24 gregl 3.5 VIEW vrev;
25 gregl 3.1 GCOORD gc[2];
26 gwlarson 3.10 FVECT cp[4], ip[4], pf, pb;
27     double af, ab, sf2, sb2, dfb2, df2, db2, penalty;
28 gregl 3.1 register int i;
29 gwlarson 3.11 /* special case */
30     if (hr <= 0 | vr <= 0)
31     return(0);
32 gregl 3.1 /* compute cell corners in image */
33     if (!hdbcoord(gc, hp, bi))
34     error(CONSISTENCY, "bad beam index in npixels");
35 gregl 3.5 hdcell(cp, hp, gc+1); /* find cell on front image */
36 gwlarson 3.10 for (i = 3; i--; ) /* compute front center */
37     pf[i] = 0.5*(cp[0][i] + cp[2][i]);
38     sf2 = 0.25*dist2(cp[0], cp[2]); /* compute half diagonal length */
39     for (i = 0; i < 4; i++) { /* compute visible quad */
40 gregl 3.5 viewloc(ip[i], vp, cp[i]);
41     if (ip[i][2] < 0.) {
42     af = 0;
43     goto getback;
44 gregl 3.4 }
45 gregl 3.5 ip[i][0] *= (double)hr; /* scale by resolution */
46     ip[i][1] *= (double)vr;
47 gregl 3.4 }
48 gregl 3.5 /* compute front area */
49     af = (ip[1][0]-ip[0][0])*(ip[2][1]-ip[0][1]) -
50     (ip[2][0]-ip[0][0])*(ip[1][1]-ip[0][1]);
51     af += (ip[2][0]-ip[3][0])*(ip[1][1]-ip[3][1]) -
52     (ip[1][0]-ip[3][0])*(ip[2][1]-ip[3][1]);
53 gwlarson 3.10 af *= af >= 0 ? 0.5 : -0.5;
54 gregl 3.5 getback:
55 schorsch 3.14 vrev = *vp; /* compute reverse view */
56 gregl 3.5 for (i = 0; i < 3; i++) {
57     vrev.vdir[i] = -vp->vdir[i];
58     vrev.vup[i] = -vp->vup[i];
59     vrev.hvec[i] = -vp->hvec[i];
60     vrev.vvec[i] = -vp->vvec[i];
61     }
62     hdcell(cp, hp, gc); /* find cell on back image */
63 gwlarson 3.10 for (i = 3; i--; ) /* compute rear center */
64     pb[i] = 0.5*(cp[0][i] + cp[2][i]);
65     sb2 = 0.25*dist2(cp[0], cp[2]); /* compute half diagonal length */
66     for (i = 0; i < 4; i++) { /* compute visible quad */
67 gregl 3.5 viewloc(ip[i], &vrev, cp[i]);
68 gwlarson 3.10 if (ip[i][2] < 0.) {
69     ab = 0;
70     goto finish;
71     }
72 gregl 3.1 ip[i][0] *= (double)hr; /* scale by resolution */
73     ip[i][1] *= (double)vr;
74     }
75 gregl 3.5 /* compute back area */
76     ab = (ip[1][0]-ip[0][0])*(ip[2][1]-ip[0][1]) -
77 gregl 3.1 (ip[2][0]-ip[0][0])*(ip[1][1]-ip[0][1]);
78 gregl 3.5 ab += (ip[2][0]-ip[3][0])*(ip[1][1]-ip[3][1]) -
79 gregl 3.1 (ip[1][0]-ip[3][0])*(ip[2][1]-ip[3][1]);
80 gwlarson 3.10 ab *= ab >= 0 ? 0.5 : -0.5;
81     finish: /* compute penalty based on dist. sightline - viewpoint */
82     df2 = dist2(vp->vp, pf);
83     db2 = dist2(vp->vp, pb);
84     dfb2 = dist2(pf, pb);
85     penalty = dfb2 + df2 - db2;
86     penalty = df2 - 0.25*penalty*penalty/dfb2;
87     if (df2 > db2) penalty /= df2 <= dfb2 ? sb2 : sb2*df2/dfb2;
88     else penalty /= db2 <= dfb2 ? sf2 : sf2*db2/dfb2;
89     if (penalty < 1.) penalty = 1.;
90     /* round off smaller non-zero area */
91     if (ab <= FTINY || (af > FTINY && af <= ab))
92     return((int)(af/penalty + 0.5));
93     return((int)(ab/penalty + 0.5));
94 gregl 3.1 }
95    
96    
97     /*
98     * The ray directions that define the pyramid in visit_cells() needn't
99     * be normalized, but they must be given in clockwise order as seen
100     * from the pyramid's apex (origin).
101 gregl 3.8 * If no cell centers fall within the domain, the closest cell is visited.
102 gregl 3.1 */
103     int
104     visit_cells(orig, pyrd, hp, vf, dp) /* visit cells within a pyramid */
105     FVECT orig, pyrd[4]; /* pyramid ray directions in clockwise order */
106 gregl 3.8 register HOLO *hp;
107 gregl 3.1 int (*vf)();
108     char *dp;
109     {
110 gregl 3.8 int ncalls = 0, n = 0;
111 gregl 3.1 int inflags = 0;
112     FVECT gp, pn[4], lo, ld;
113     double po[4], lbeg, lend, d, t;
114 gregl 3.8 GCOORD gc, gc2[2];
115 gregl 3.1 register int i;
116     /* figure out whose side we're on */
117     hdgrid(gp, hp, orig);
118     for (i = 0; i < 3; i++) {
119     inflags |= (gp[i] > FTINY) << (i<<1);
120     inflags |= (gp[i] < hp->grid[i]-FTINY) << (i<<1 | 1);
121     }
122     /* compute pyramid planes */
123     for (i = 0; i < 4; i++) {
124     fcross(pn[i], pyrd[i], pyrd[(i+1)&03]);
125     po[i] = DOT(pn[i], orig);
126     }
127     /* traverse each wall */
128     for (gc.w = 0; gc.w < 6; gc.w++) {
129     if (!(inflags & 1<<gc.w)) /* origin on wrong side */
130     continue;
131     /* scanline algorithm */
132 gregl 3.9 for (gc.i[1] = hp->grid[hdwg1[gc.w]]; gc.i[1]--; ) {
133 gregl 3.1 /* compute scanline */
134     gp[gc.w>>1] = gc.w&1 ? hp->grid[gc.w>>1] : 0;
135 gregl 3.9 gp[hdwg0[gc.w]] = 0;
136     gp[hdwg1[gc.w]] = gc.i[1] + 0.5;
137 gregl 3.1 hdworld(lo, hp, gp);
138 gregl 3.9 gp[hdwg0[gc.w]] = 1;
139 gregl 3.1 hdworld(ld, hp, gp);
140 gregl 3.2 ld[0] -= lo[0]; ld[1] -= lo[1]; ld[2] -= lo[2];
141 gregl 3.1 /* find scanline limits */
142 gregl 3.9 lbeg = 0; lend = hp->grid[hdwg0[gc.w]];
143 gregl 3.1 for (i = 0; i < 4; i++) {
144     t = DOT(pn[i], lo) - po[i];
145     d = -DOT(pn[i], ld);
146 gregl 3.2 if (d > FTINY) { /* <- plane */
147 gregl 3.1 if ((t /= d) < lend)
148     lend = t;
149 gregl 3.2 } else if (d < -FTINY) { /* plane -> */
150 gregl 3.1 if ((t /= d) > lbeg)
151     lbeg = t;
152 gregl 3.3 } else if (t < 0) { /* outside */
153     lend = -1;
154     break;
155     }
156 gregl 3.1 }
157 gregl 3.3 if (lbeg >= lend)
158     continue;
159 gregl 3.1 i = lend + .5; /* visit cells on this scan */
160 gregl 3.8 for (gc.i[0] = lbeg + .5; gc.i[0] < i; gc.i[0]++) {
161 gregl 3.1 n += (*vf)(&gc, dp);
162 gregl 3.8 ncalls++;
163     }
164 gregl 3.1 }
165     }
166 gregl 3.8 if (ncalls) /* got one at least */
167     return(n);
168     /* else find closest cell */
169     VSUM(ld, pyrd[0], pyrd[1], 1.);
170     VSUM(ld, ld, pyrd[2], 1.);
171     VSUM(ld, ld, pyrd[3], 1.);
172     #if 0
173     if (normalize(ld) == 0.0) /* technically not necessary */
174     return(0);
175     #endif
176     d = hdinter(gc2, NULL, &t, hp, orig, ld);
177     if (d >= FHUGE || t <= 0.)
178     return(0);
179     return((*vf)(gc2+1, dp)); /* visit it */
180 gregl 3.1 }
181    
182    
183 gregl 3.4 sect_behind(hp, vp) /* check if section is "behind" viewpoint */
184     register HOLO *hp;
185     register VIEW *vp;
186     {
187     FVECT hcent;
188     /* compute holodeck section center */
189     VSUM(hcent, hp->orig, hp->xv[0], 0.5);
190     VSUM(hcent, hcent, hp->xv[1], 0.5);
191     VSUM(hcent, hcent, hp->xv[2], 0.5);
192     /* behind if center is behind */
193     return(DOT(vp->vdir,hcent) < DOT(vp->vdir,vp->vp));
194     }
195    
196    
197     viewpyramid(org, dir, hp, vp) /* compute view pyramid */
198     FVECT org, dir[4];
199     HOLO *hp;
200     VIEW *vp;
201     {
202     register int i;
203     /* check view type */
204     if (vp->type == VT_PAR)
205     return(0);
206     /* in front or behind? */
207     if (!sect_behind(hp, vp)) {
208     if (viewray(org, dir[0], vp, 0., 0.) < -FTINY)
209     return(0);
210     if (viewray(org, dir[1], vp, 0., 1.) < -FTINY)
211     return(0);
212     if (viewray(org, dir[2], vp, 1., 1.) < -FTINY)
213     return(0);
214     if (viewray(org, dir[3], vp, 1., 0.) < -FTINY)
215     return(0);
216     return(1);
217     } /* reverse pyramid */
218     if (viewray(org, dir[3], vp, 0., 0.) < -FTINY)
219     return(0);
220     if (viewray(org, dir[2], vp, 0., 1.) < -FTINY)
221     return(0);
222     if (viewray(org, dir[1], vp, 1., 1.) < -FTINY)
223     return(0);
224     if (viewray(org, dir[0], vp, 1., 0.) < -FTINY)
225     return(0);
226     for (i = 0; i < 3; i++) {
227     dir[0][i] = -dir[0][i];
228     dir[1][i] = -dir[1][i];
229     dir[2][i] = -dir[2][i];
230     dir[3][i] = -dir[3][i];
231     }
232     return(-1);
233     }
234    
235    
236 gregl 3.1 int
237     addcell(gcp, cl) /* add a cell to a list */
238     GCOORD *gcp;
239 gregl 3.4 register struct cellist *cl;
240 gregl 3.1 {
241 schorsch 3.14 *(cl->cl+cl->n) = *gcp;
242 gregl 3.4 cl->n++;
243 gregl 3.1 return(1);
244     }
245    
246    
247     int
248     cellcmp(gcp1, gcp2) /* visit_cells() cell ordering */
249     register GCOORD *gcp1, *gcp2;
250     {
251     register int c;
252    
253     if ((c = gcp1->w - gcp2->w))
254     return(c);
255     if ((c = gcp2->i[1] - gcp1->i[1])) /* wg1 is reverse-ordered */
256     return(c);
257     return(gcp1->i[0] - gcp2->i[0]);
258     }
259    
260    
261 gregl 3.4 GCOORD *
262     getviewcells(np, hp, vp) /* get ordered cell list for section view */
263     int *np; /* returned number of cells (negative if reversed) */
264 gregl 3.1 register HOLO *hp;
265     VIEW *vp;
266     {
267     FVECT org, dir[4];
268 gregl 3.4 int orient;
269     struct cellist cl;
270 gregl 3.1 /* compute view pyramid */
271 gregl 3.4 *np = 0;
272     orient = viewpyramid(org, dir, hp, vp);
273     if (!orient)
274     return(NULL);
275 gregl 3.1 /* allocate enough list space */
276 gregl 3.4 cl.n = 2*( hp->grid[0]*hp->grid[1] +
277     hp->grid[0]*hp->grid[2] +
278     hp->grid[1]*hp->grid[2] );
279     cl.cl = (GCOORD *)malloc(cl.n*sizeof(GCOORD));
280     if (cl.cl == NULL)
281 gregl 3.1 goto memerr;
282 gregl 3.4 cl.n = 0; /* add cells within pyramid */
283 gwlarson 3.11 visit_cells(org, dir, hp, addcell, (char *)&cl);
284 gregl 3.4 if (!cl.n) {
285 greg 3.12 free((void *)cl.cl);
286 gregl 3.1 return(NULL);
287     }
288 gregl 3.4 *np = cl.n * orient;
289 gregl 3.1 #if 0
290 gregl 3.2 /* We're just going to free this memory in a moment, and list is
291     * sorted automatically by visit_cells(), so we don't need this.
292     */
293 gregl 3.4 /* optimize memory use */
294 greg 3.13 cl.cl = (GCOORD *)realloc((void *)cl.cl, cl.n*sizeof(GCOORD));
295 gregl 3.4 if (cl.cl == NULL)
296     goto memerr;
297 gregl 3.1 /* sort the list */
298 gregl 3.4 qsort((char *)cl.cl, cl.n, sizeof(GCOORD), cellcmp);
299 gregl 3.1 #endif
300 gregl 3.4 return(cl.cl);
301 gregl 3.1 memerr:
302     error(SYSTEM, "out of memory in getviewcells");
303     }
304 gregl 3.6
305    
306 schorsch 3.15 extern void
307     gridlines( /* run through holodeck section grid lines */
308     void (*f)(FVECT wp[2])
309     )
310 gregl 3.6 {
311     register int hd, w, i;
312     int g0, g1;
313 gregl 3.7 FVECT wp[2], mov;
314 gregl 3.6 double d;
315     /* do each wall on each section */
316     for (hd = 0; hdlist[hd] != NULL; hd++)
317     for (w = 0; w < 6; w++) {
318 gregl 3.9 g0 = hdwg0[w];
319     g1 = hdwg1[w];
320 gregl 3.7 d = 1.0/hdlist[hd]->grid[g0];
321     mov[0] = d * hdlist[hd]->xv[g0][0];
322     mov[1] = d * hdlist[hd]->xv[g0][1];
323     mov[2] = d * hdlist[hd]->xv[g0][2];
324     if (w & 1) {
325 gregl 3.6 VSUM(wp[0], hdlist[hd]->orig,
326     hdlist[hd]->xv[w>>1], 1.);
327 gregl 3.7 VSUM(wp[0], wp[0], mov, 1.);
328     } else
329     VCOPY(wp[0], hdlist[hd]->orig);
330     VSUM(wp[1], wp[0], hdlist[hd]->xv[g1], 1.);
331     for (i = hdlist[hd]->grid[g0]; ; ) { /* g0 lines */
332 gregl 3.6 (*f)(wp);
333 gregl 3.7 if (!--i) break;
334     wp[0][0] += mov[0]; wp[0][1] += mov[1];
335     wp[0][2] += mov[2]; wp[1][0] += mov[0];
336     wp[1][1] += mov[1]; wp[1][2] += mov[2];
337 gregl 3.6 }
338 gregl 3.7 d = 1.0/hdlist[hd]->grid[g1];
339     mov[0] = d * hdlist[hd]->xv[g1][0];
340     mov[1] = d * hdlist[hd]->xv[g1][1];
341     mov[2] = d * hdlist[hd]->xv[g1][2];
342     if (w & 1)
343 gregl 3.6 VSUM(wp[0], hdlist[hd]->orig,
344     hdlist[hd]->xv[w>>1], 1.);
345 gregl 3.7 else
346     VSUM(wp[0], hdlist[hd]->orig, mov, 1.);
347     VSUM(wp[1], wp[0], hdlist[hd]->xv[g0], 1.);
348     for (i = hdlist[hd]->grid[g1]; ; ) { /* g1 lines */
349 gregl 3.6 (*f)(wp);
350 gregl 3.7 if (!--i) break;
351     wp[0][0] += mov[0]; wp[0][1] += mov[1];
352     wp[0][2] += mov[2]; wp[1][0] += mov[0];
353     wp[1][1] += mov[1]; wp[1][2] += mov[2];
354 gregl 3.6 }
355     }
356     }