ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/hd/rhpict2.c
Revision: 3.5
Committed: Tue Mar 9 11:41:29 1999 UTC (25 years, 1 month ago) by gwlarson
Content type: text/plain
Branch: MAIN
Changes since 3.4: +9 -6 lines
Log Message:
bug fix in occlusion culling

File Contents

# Content
1 /* Copyright (c) 1999 Silicon Graphics, Inc. */
2
3 #ifndef lint
4 static char SCCSid[] = "$SunId$ SGI";
5 #endif
6
7 /*
8 * Rendering routines for rhpict.
9 */
10
11 #include "holo.h"
12 #include "view.h"
13
14 #ifndef DEPS
15 #define DEPS 0.02 /* depth epsilon */
16 #endif
17 #ifndef PEPS
18 #define PEPS 0.04 /* pixel value epsilon */
19 #endif
20 #ifndef MAXRAD
21 #define MAXRAD 64 /* maximum kernel radius */
22 #endif
23 #ifndef NNEIGH
24 #define NNEIGH 7 /* find this many neighbors */
25 #endif
26
27 #define NINF 16382
28
29 #define MAXRAD2 (MAXRAD*MAXRAD+1)
30
31 #define G0NORM 0.286 /* ground zero normalization (1/x integral) */
32
33 #ifndef FL4OP
34 #define FL4OP(f,i,op) ((f)[(i)>>5] op (1L<<((i)&0x1f)))
35 #define CHK4(f,i) FL4OP(f,i,&)
36 #define SET4(f,i) FL4OP(f,i,|=)
37 #define CLR4(f,i) FL4OP(f,i,&=~)
38 #define TGL4(f,i) FL4OP(f,i,^=)
39 #define FL4NELS(n) (((n)+0x1f)>>5)
40 #define CLR4ALL(f,n) bzero((char *)(f),FL4NELS(n)*sizeof(int4))
41 #endif
42
43 static int4 *pixFlags; /* pixel occupancy flags */
44 static float pixWeight[MAXRAD2]; /* pixel weighting function */
45 static short isqrttab[MAXRAD2]; /* integer square root table */
46
47 #define isqrt(i2) ((int)isqrttab[(int)(i2)])
48
49 extern VIEW myview; /* current output view */
50 extern COLOR *mypixel; /* pixels being rendered */
51 extern float *myweight; /* weights (used to compute final pixels) */
52 extern float *mydepth; /* depth values (visibility culling) */
53 extern int hres, vres; /* current horizontal and vertical res. */
54
55
56 pixBeam(bp, hb) /* render a particular beam */
57 BEAM *bp;
58 register HDBEAMI *hb;
59 {
60 GCOORD gc[2];
61 register RAYVAL *rv;
62 FVECT rorg, rdir, wp, ip;
63 double d, prox;
64 COLOR col;
65 int n;
66 register int4 p;
67
68 if (!hdbcoord(gc, hb->h, hb->b))
69 error(CONSISTENCY, "bad beam in render_beam");
70 for (n = bp->nrm, rv = hdbray(bp); n--; rv++) {
71 /* reproject each sample */
72 hdray(rorg, rdir, hb->h, gc, rv->r);
73 if (rv->d < DCINF) {
74 d = hddepth(hb->h, rv->d);
75 VSUM(wp, rorg, rdir, d);
76 VSUB(ip, wp, myview.vp);
77 d = DOT(ip,rdir);
78 prox = d*d/DOT(ip,ip); /* cos(diff_angle)^32 */
79 prox *= prox; prox *= prox; prox *= prox; prox *= prox;
80 } else {
81 if (myview.type == VT_PAR || myview.vaft > FTINY)
82 continue; /* inf. off view */
83 VSUM(wp, myview.vp, rdir, FHUGE);
84 prox = 1.;
85 }
86 viewloc(ip, &myview, wp); /* frustum clipping */
87 if (ip[2] < 0.)
88 continue;
89 if (ip[0] < 0. || ip[0] >= 1.)
90 continue;
91 if (ip[1] < 0. || ip[1] >= 1.)
92 continue;
93 if (myview.vaft > FTINY && ip[2] > myview.vaft - myview.vfore)
94 continue; /* not exact for VT_PER */
95 p = (int)(ip[1]*vres)*hres + (int)(ip[0]*hres);
96 if (mydepth[p] > FTINY) { /* check depth */
97 if (ip[2] > mydepth[p]*(1.+DEPS))
98 continue;
99 if (ip[2] < mydepth[p]*(1.-DEPS)) {
100 setcolor(mypixel[p], 0., 0., 0.);
101 myweight[p] = 0.;
102 }
103 }
104 colr_color(col, rv->v);
105 scalecolor(col, prox);
106 addcolor(mypixel[p], col);
107 myweight[p] += prox;
108 mydepth[p] = ip[2];
109 }
110 }
111
112
113 int
114 kill_occl(h, v, nl, n) /* check for occlusion errors */
115 int h, v;
116 register short nl[NNEIGH][2];
117 int n;
118 {
119 short forequad[2][2];
120 int d;
121 register int4 i, p;
122
123 if (n <= 0) {
124 #ifdef DEBUG
125 error(WARNING, "neighborless sample in kill_occl");
126 #endif
127 return(1);
128 }
129 p = v*hres + h;
130 forequad[0][0] = forequad[0][1] = forequad[1][0] = forequad[1][1] = 0;
131 for (i = n; i--; ) {
132 d = (h-nl[i][0])*(h-nl[i][0]) + (v-nl[i][1])*(v-nl[i][1]);
133 d = isqrt(d);
134 if (mydepth[nl[i][1]*hres+nl[i][0]]*(1.+DEPS*d) < mydepth[p])
135 forequad[nl[i][0]<h][nl[i][1]<v] = 1;
136 }
137 if (forequad[0][0]+forequad[0][1]+forequad[1][0]+forequad[1][1] > 2) {
138 setcolor(mypixel[p], 0., 0., 0.);
139 myweight[p] = 0.; /* occupancy reset afterwards */
140 }
141 return(1);
142 }
143
144
145 int
146 grow_samp(h, v, nl, n) /* grow sample point appropriately */
147 int h, v;
148 register short nl[NNEIGH][2];
149 int n;
150 {
151 int dis[NNEIGH], ndis;
152 COLOR mykern[MAXRAD2];
153 int4 maxr2;
154 double d;
155 register int4 p, r2;
156 int i, r, maxr, h2, v2;
157
158 if (n <= 0)
159 return(1);
160 p = v*hres + h; /* build kernel values */
161 maxr2 = (h-nl[n-1][0])*(h-nl[n-1][0]) + (v-nl[n-1][1])*(v-nl[n-1][1]);
162 DCHECK(maxr2>=MAXRAD2, CONSISTENCY, "out of range neighbor");
163 maxr = isqrt(maxr2);
164 for (v2 = 1; v2 <= maxr; v2++)
165 for (h2 = 0; h2 <= v2; h2++) {
166 r2 = h2*h2 + v2*v2;
167 if (r2 > maxr2) break;
168 copycolor(mykern[r2], mypixel[p]);
169 scalecolor(mykern[r2], pixWeight[r2]);
170 }
171 ndis = 0; /* find discontinuities */
172 for (i = n; i--; ) {
173 r2 = (h-nl[i][0])*(h-nl[i][0]) + (v-nl[i][1])*(v-nl[i][1]);
174 r = isqrt(r2);
175 d = mydepth[nl[i][1]*hres+nl[i][0]] / mydepth[p];
176 d = d>=1. ? d-1. : 1.-d;
177 if (d > r*DEPS || bigdiff(mypixel[p],
178 mypixel[nl[i][1]*hres+nl[i][0]], r*PEPS))
179 dis[ndis++] = i;
180 }
181 /* stamp out that kernel */
182 for (v2 = v-maxr; v2 <= v+maxr; v2++) {
183 if (v2 < 0) v2 = 0;
184 else if (v2 >= vres) break;
185 for (h2 = h-maxr; h2 <= h+maxr; h2++) {
186 if (h2 < 0) h2 = 0;
187 else if (h2 >= hres) break;
188 r2 = (h2-h)*(h2-h) + (v2-v)*(v2-v);
189 if (r2 > maxr2) continue;
190 if (CHK4(pixFlags, v2*hres+h2))
191 continue; /* occupied */
192 for (i = ndis; i--; ) {
193 r = (h2-nl[dis[i]][0])*(h2-nl[dis[i]][0]) +
194 (v2-nl[dis[i]][1])*(v2-nl[dis[i]][1]);
195 if (r < r2) break;
196 }
197 if (i >= 0) continue; /* outside edge */
198 addcolor(mypixel[v2*hres+h2], mykern[r2]);
199 myweight[v2*hres+h2] += pixWeight[r2] * myweight[p];
200 }
201 }
202 return(1);
203 }
204
205
206 pixFlush() /* done with beams -- flush pixel values */
207 {
208 if (pixWeight[0] <= FTINY)
209 init_wfunc(); /* initialize weighting function */
210 reset_flags(); /* set occupancy flags */
211 meet_neighbors(kill_occl); /* eliminate occlusion errors */
212 reset_flags(); /* reset occupancy flags */
213 meet_neighbors(grow_samp); /* grow valid samples over image */
214 free((char *)pixFlags); /* free pixel flags */
215 pixFlags = NULL;
216 }
217
218
219 reset_flags() /* allocate/set/reset occupancy flags */
220 {
221 register int4 p;
222
223 if (pixFlags == NULL) {
224 pixFlags = (int4 *)calloc(FL4NELS(hres*vres), sizeof(int4));
225 CHECK(pixFlags==NULL, SYSTEM, "out of memory in reset_flags");
226 } else
227 CLR4ALL(pixFlags, hres*vres);
228 for (p = hres*vres; p--; )
229 if (myweight[p] > FTINY)
230 SET4(pixFlags, p);
231 }
232
233
234 init_wfunc() /* initialize weighting function */
235 {
236 register int i, j;
237 register int4 r2;
238 register double d;
239
240 for (i = 1; i <= MAXRAD; i++)
241 for (j = 0; j <= i; j++) {
242 r2 = i*i + j*j;
243 if (r2 >= MAXRAD2) break;
244 d = sqrt((double)r2);
245 pixWeight[r2] = G0NORM/d;
246 isqrttab[r2] = d + 0.99;
247 }
248 pixWeight[0] = 1.;
249 isqrttab[0] = 0;
250 }
251
252
253 int
254 findneigh(nl, h, v, rnl) /* find NNEIGH neighbors for pixel */
255 short nl[NNEIGH][2];
256 int h, v;
257 register short (*rnl)[NNEIGH];
258 {
259 int nn = 0;
260 int4 d, nd[NNEIGH];
261 int n, hoff;
262 register int h2, n2;
263
264 nd[NNEIGH-1] = MAXRAD2;
265 for (hoff = 1; hoff < hres; hoff = (hoff<0) - hoff) {
266 h2 = h + hoff;
267 if (h2 < 0 | h2 >= hres)
268 continue;
269 if ((h2-h)*(h2-h) >= nd[NNEIGH-1])
270 break;
271 for (n = 0; n < NNEIGH && rnl[h2][n] < NINF; n++) {
272 d = (h2-h)*(h2-h) + (v-rnl[h2][n])*(v-rnl[h2][n]);
273 if (d >= nd[NNEIGH-1])
274 continue;
275 if (nn < NNEIGH) /* insert neighbor */
276 nn++;
277 for (n2 = nn; n2--; ) {
278 if (!n2 || d >= nd[n2-1]) {
279 nd[n2] = d;
280 nl[n2][0] = h2;
281 nl[n2][1] = rnl[h2][n];
282 break;
283 }
284 nd[n2] = nd[n2-1];
285 nl[n2][0] = nl[n2-1][0];
286 nl[n2][1] = nl[n2-1][1];
287 }
288 }
289 }
290 return(nn);
291 }
292
293
294 meet_neighbors(nf) /* run through samples and their neighbors */
295 int (*nf)();
296 {
297 short ln[NNEIGH][2];
298 int h, v, n, v2;
299 register short (*rnl)[NNEIGH];
300 /* initialize bottom row list */
301 rnl = (short (*)[NNEIGH])malloc(NNEIGH*sizeof(short)*hres);
302 CHECK(rnl==NULL, SYSTEM, "out of memory in meet_neighbors");
303 for (h = 0; h < hres; h++) {
304 for (n = v = 0; v < vres; v++)
305 if (CHK4(pixFlags, v*hres+h)) {
306 rnl[h][n++] = v;
307 if (n >= NNEIGH)
308 break;
309 }
310 while (n < NNEIGH)
311 rnl[h][n++] = NINF;
312 }
313 v = 0; /* do each row */
314 for ( ; ; ) {
315 for (h = 0; h < hres; h++) {
316 if (!CHK4(pixFlags, v*hres+h))
317 continue; /* no one home */
318 n = findneigh(ln, h, v, rnl);
319 (*nf)(h, v, ln, n); /* call on neighbors */
320 }
321 if (++v >= vres) /* reinitialize row list */
322 break;
323 for (h = 0; h < hres; h++)
324 for (v2 = rnl[h][NNEIGH-1]+1; v2 < vres; v2++) {
325 if (v2 - v > v - rnl[h][0])
326 break; /* not close enough */
327 if (CHK4(pixFlags, v2*hres+h)) {
328 for (n = 0; n < NNEIGH-1; n++)
329 rnl[h][n] = rnl[h][n+1];
330 rnl[h][NNEIGH-1] = v2;
331 }
332 }
333 }
334 free((char *)rnl); /* free row list */
335 }