ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/hd/rhpict2.c
Revision: 3.15
Committed: Thu Jan 1 11:21:55 2004 UTC (20 years, 3 months ago) by schorsch
Content type: text/plain
Branch: MAIN
CVS Tags: rad4R2P2, rad5R0, rad5R1, rad3R7P2, rad3R7P1, rad4R2, rad4R1, rad4R0, rad3R6, rad3R6P1, rad3R8, rad3R9, rad4R2P1
Changes since 3.14: +69 -38 lines
Log Message:
Ansification and prototypes.

File Contents

# User Rev Content
1 gwlarson 3.1 #ifndef lint
2 schorsch 3.15 static const char RCSid[] = "$Id: rhpict2.c,v 3.14 2003/07/27 22:12:02 schorsch Exp $";
3 gwlarson 3.1 #endif
4     /*
5     * Rendering routines for rhpict.
6     */
7    
8 schorsch 3.12 #include <string.h>
9    
10 gwlarson 3.1 #include "holo.h"
11     #include "view.h"
12 gwlarson 3.6 #include "random.h"
13 gwlarson 3.1
14 gwlarson 3.2 #ifndef DEPS
15 gwlarson 3.3 #define DEPS 0.02 /* depth epsilon */
16 gwlarson 3.2 #endif
17 gwlarson 3.3 #ifndef PEPS
18     #define PEPS 0.04 /* pixel value epsilon */
19     #endif
20 gwlarson 3.2 #ifndef MAXRAD
21     #define MAXRAD 64 /* maximum kernel radius */
22     #endif
23     #ifndef NNEIGH
24 gwlarson 3.9 #define NNEIGH 5 /* find this many neighbors */
25 gwlarson 3.2 #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 schorsch 3.12 #define CLR4ALL(f,n) memset((char *)(f),'\0',FL4NELS(n)*sizeof(int32))
41 gwlarson 3.2 #endif
42    
43 greg 3.11 static int32 *pixFlags; /* pixel occupancy flags */
44 gwlarson 3.2 static float pixWeight[MAXRAD2]; /* pixel weighting function */
45 gwlarson 3.3 static short isqrttab[MAXRAD2]; /* integer square root table */
46 gwlarson 3.2
47 gwlarson 3.9 #define isqrt(i2) (isqrttab[i2])
48 gwlarson 3.3
49 gwlarson 3.1 extern VIEW myview; /* current output view */
50     extern COLOR *mypixel; /* pixels being rendered */
51     extern float *myweight; /* weights (used to compute final pixels) */
52 gwlarson 3.2 extern float *mydepth; /* depth values (visibility culling) */
53 gwlarson 3.1 extern int hres, vres; /* current horizontal and vertical res. */
54    
55 schorsch 3.15 extern void pixFinish(double ransamp);
56     extern void pixBeam(BEAM *bp, HDBEAMI *hb);
57 gwlarson 3.1
58 schorsch 3.15 typedef int (sampfunc)(int h, int v, short nl[NNEIGH][2], int nd[NNEIGH],
59     int n, double *rf);
60     static sampfunc kill_occl;
61     static sampfunc smooth_samp;
62     static sampfunc random_samp;
63     static void meet_neighbors(int occ, sampfunc *nf, double *dp);
64     static void reset_flags(void);
65     static void init_wfunc(void);
66     static int findneigh(short nl[NNEIGH][2], int nd[NNEIGH], int h, int v,
67     register short (*rnl)[NNEIGH]);
68    
69    
70     extern void
71     pixBeam( /* render a particular beam */
72     BEAM *bp,
73     register HDBEAMI *hb
74     )
75 gwlarson 3.1 {
76     GCOORD gc[2];
77     register RAYVAL *rv;
78     FVECT rorg, rdir, wp, ip;
79     double d, prox;
80     COLOR col;
81 gwlarson 3.2 int n;
82 greg 3.11 register int32 p;
83 gwlarson 3.1
84     if (!hdbcoord(gc, hb->h, hb->b))
85     error(CONSISTENCY, "bad beam in render_beam");
86     for (n = bp->nrm, rv = hdbray(bp); n--; rv++) {
87     /* reproject each sample */
88     hdray(rorg, rdir, hb->h, gc, rv->r);
89     if (rv->d < DCINF) {
90     d = hddepth(hb->h, rv->d);
91     VSUM(wp, rorg, rdir, d);
92     VSUB(ip, wp, myview.vp);
93     d = DOT(ip,rdir);
94 gwlarson 3.2 prox = d*d/DOT(ip,ip); /* cos(diff_angle)^32 */
95     prox *= prox; prox *= prox; prox *= prox; prox *= prox;
96 gwlarson 3.1 } else {
97     if (myview.type == VT_PAR || myview.vaft > FTINY)
98     continue; /* inf. off view */
99     VSUM(wp, myview.vp, rdir, FHUGE);
100 gwlarson 3.2 prox = 1.;
101 gwlarson 3.1 }
102     viewloc(ip, &myview, wp); /* frustum clipping */
103     if (ip[2] < 0.)
104     continue;
105     if (ip[0] < 0. || ip[0] >= 1.)
106     continue;
107     if (ip[1] < 0. || ip[1] >= 1.)
108     continue;
109     if (myview.vaft > FTINY && ip[2] > myview.vaft - myview.vfore)
110 gwlarson 3.2 continue; /* not exact for VT_PER */
111 gwlarson 3.1 p = (int)(ip[1]*vres)*hres + (int)(ip[0]*hres);
112 gwlarson 3.2 if (mydepth[p] > FTINY) { /* check depth */
113     if (ip[2] > mydepth[p]*(1.+DEPS))
114     continue;
115     if (ip[2] < mydepth[p]*(1.-DEPS)) {
116     setcolor(mypixel[p], 0., 0., 0.);
117     myweight[p] = 0.;
118     }
119     }
120 gwlarson 3.1 colr_color(col, rv->v);
121 gwlarson 3.2 scalecolor(col, prox);
122 gwlarson 3.1 addcolor(mypixel[p], col);
123 gwlarson 3.2 myweight[p] += prox;
124     mydepth[p] = ip[2];
125 gwlarson 3.1 }
126 gwlarson 3.2 }
127    
128    
129 schorsch 3.15 static int
130     kill_occl( /* check for occlusion errors */
131     int h,
132     int v,
133     register short nl[NNEIGH][2],
134     int nd[NNEIGH],
135     int n,
136     double *rf
137     )
138 gwlarson 3.2 {
139     short forequad[2][2];
140     int d;
141 gwlarson 3.9 register int i;
142 greg 3.11 register int32 p;
143 gwlarson 3.2
144 gwlarson 3.5 if (n <= 0) {
145     #ifdef DEBUG
146     error(WARNING, "neighborless sample in kill_occl");
147     #endif
148 gwlarson 3.2 return(1);
149 gwlarson 3.5 }
150 gwlarson 3.3 p = v*hres + h;
151 gwlarson 3.2 forequad[0][0] = forequad[0][1] = forequad[1][0] = forequad[1][1] = 0;
152     for (i = n; i--; ) {
153 gwlarson 3.9 d = isqrt(nd[i]);
154 gwlarson 3.3 if (mydepth[nl[i][1]*hres+nl[i][0]]*(1.+DEPS*d) < mydepth[p])
155 gwlarson 3.2 forequad[nl[i][0]<h][nl[i][1]<v] = 1;
156     }
157 gwlarson 3.3 if (forequad[0][0]+forequad[0][1]+forequad[1][0]+forequad[1][1] > 2) {
158     setcolor(mypixel[p], 0., 0., 0.);
159     myweight[p] = 0.; /* occupancy reset afterwards */
160 gwlarson 3.2 }
161     return(1);
162     }
163    
164    
165 schorsch 3.15 static int
166     smooth_samp( /* grow sample point smoothly */
167     int h,
168     int v,
169     register short nl[NNEIGH][2],
170     int nd[NNEIGH],
171     int n,
172     double *rf
173     )
174 gwlarson 3.2 {
175 gwlarson 3.3 int dis[NNEIGH], ndis;
176 gwlarson 3.2 COLOR mykern[MAXRAD2];
177 gwlarson 3.9 int maxr2;
178 gwlarson 3.5 double d;
179 greg 3.11 register int32 p;
180 gwlarson 3.9 register int r2;
181 gwlarson 3.3 int i, r, maxr, h2, v2;
182 gwlarson 3.2
183     if (n <= 0)
184     return(1);
185     p = v*hres + h; /* build kernel values */
186 gwlarson 3.9 maxr2 = nd[n-1];
187 gwlarson 3.2 DCHECK(maxr2>=MAXRAD2, CONSISTENCY, "out of range neighbor");
188 gwlarson 3.3 maxr = isqrt(maxr2);
189     for (v2 = 1; v2 <= maxr; v2++)
190     for (h2 = 0; h2 <= v2; h2++) {
191     r2 = h2*h2 + v2*v2;
192     if (r2 > maxr2) break;
193     copycolor(mykern[r2], mypixel[p]);
194     scalecolor(mykern[r2], pixWeight[r2]);
195     }
196     ndis = 0; /* find discontinuities */
197     for (i = n; i--; ) {
198 gwlarson 3.9 r = isqrt(nd[i]);
199 gwlarson 3.3 d = mydepth[nl[i][1]*hres+nl[i][0]] / mydepth[p];
200     d = d>=1. ? d-1. : 1.-d;
201     if (d > r*DEPS || bigdiff(mypixel[p],
202     mypixel[nl[i][1]*hres+nl[i][0]], r*PEPS))
203     dis[ndis++] = i;
204 gwlarson 3.2 }
205 gwlarson 3.3 /* stamp out that kernel */
206 gwlarson 3.2 for (v2 = v-maxr; v2 <= v+maxr; v2++) {
207 gwlarson 3.3 if (v2 < 0) v2 = 0;
208     else if (v2 >= vres) break;
209 gwlarson 3.2 for (h2 = h-maxr; h2 <= h+maxr; h2++) {
210 gwlarson 3.3 if (h2 < 0) h2 = 0;
211     else if (h2 >= hres) break;
212     r2 = (h2-h)*(h2-h) + (v2-v)*(v2-v);
213 gwlarson 3.2 if (r2 > maxr2) continue;
214     if (CHK4(pixFlags, v2*hres+h2))
215     continue; /* occupied */
216 gwlarson 3.3 for (i = ndis; i--; ) {
217     r = (h2-nl[dis[i]][0])*(h2-nl[dis[i]][0]) +
218     (v2-nl[dis[i]][1])*(v2-nl[dis[i]][1]);
219     if (r < r2) break;
220     }
221     if (i >= 0) continue; /* outside edge */
222 gwlarson 3.2 addcolor(mypixel[v2*hres+h2], mykern[r2]);
223 gwlarson 3.5 myweight[v2*hres+h2] += pixWeight[r2] * myweight[p];
224 gwlarson 3.2 }
225     }
226     return(1);
227     }
228    
229    
230 schorsch 3.15 static int
231     random_samp( /* gather samples randomly */
232     int h,
233     int v,
234     register short nl[NNEIGH][2],
235     int nd[NNEIGH],
236     int n,
237     double *rf
238     )
239 gwlarson 3.2 {
240 gwlarson 3.9 float rnt[NNEIGH];
241     double rvar;
242 greg 3.11 register int32 p, pn;
243 gwlarson 3.9 register int ni;
244 gwlarson 3.6
245     if (n <= 0)
246     return(1);
247 gwlarson 3.9 p = v*hres + h;
248     if (*rf <= FTINY) /* straight Voronoi regions */
249     ni = 0;
250     else { /* weighted choice */
251     DCHECK(nd[n-1]>=MAXRAD2, CONSISTENCY, "out of range neighbor");
252     rnt[0] = pixWeight[nd[0]];
253     for (ni = 1; ni < n; ni++)
254     rnt[ni] = rnt[ni-1] + pixWeight[nd[ni]];
255     rvar = rnt[n-1]*pow(frandom(), 1. / *rf);
256     for (ni = 0; rvar > rnt[ni]+FTINY; ni++)
257     ;
258 gwlarson 3.6 }
259 gwlarson 3.9 pn = nl[ni][1]*hres + nl[ni][0];
260     addcolor(mypixel[p], mypixel[pn]);
261     myweight[p] += myweight[pn];
262 gwlarson 3.6 return(1);
263     }
264    
265    
266 schorsch 3.15 extern void
267     pixFinish( /* done with beams -- compute pixel values */
268     double ransamp
269     )
270 gwlarson 3.6 {
271 gwlarson 3.5 if (pixWeight[0] <= FTINY)
272     init_wfunc(); /* initialize weighting function */
273 gwlarson 3.2 reset_flags(); /* set occupancy flags */
274 gwlarson 3.9 meet_neighbors(1,kill_occl,NULL); /* identify occlusion errors */
275 gwlarson 3.2 reset_flags(); /* reset occupancy flags */
276 gwlarson 3.7 if (ransamp >= 0.) /* spread samples over image */
277 gwlarson 3.9 meet_neighbors(0,random_samp,&ransamp);
278 gwlarson 3.6 else
279 gwlarson 3.9 meet_neighbors(1,smooth_samp,NULL);
280 greg 3.10 free((void *)pixFlags); /* free pixel flags */
281 gwlarson 3.2 pixFlags = NULL;
282     }
283    
284    
285 schorsch 3.15 static void
286     reset_flags(void) /* allocate/set/reset occupancy flags */
287 gwlarson 3.2 {
288 greg 3.11 register int32 p;
289 gwlarson 3.2
290     if (pixFlags == NULL) {
291 greg 3.11 pixFlags = (int32 *)calloc(FL4NELS(hres*vres), sizeof(int32));
292 gwlarson 3.2 CHECK(pixFlags==NULL, SYSTEM, "out of memory in reset_flags");
293     } else
294     CLR4ALL(pixFlags, hres*vres);
295     for (p = hres*vres; p--; )
296     if (myweight[p] > FTINY)
297     SET4(pixFlags, p);
298 gwlarson 3.4 }
299    
300    
301 schorsch 3.15 static void
302     init_wfunc(void) /* initialize weighting function */
303 gwlarson 3.4 {
304 gwlarson 3.9 register int r2;
305 gwlarson 3.4 register double d;
306    
307 gwlarson 3.6 for (r2 = MAXRAD2; --r2; ) {
308     d = sqrt((double)r2);
309     pixWeight[r2] = G0NORM/d;
310     isqrttab[r2] = d + 0.99;
311     }
312 gwlarson 3.4 pixWeight[0] = 1.;
313     isqrttab[0] = 0;
314 gwlarson 3.2 }
315    
316    
317 schorsch 3.15 static int
318     findneigh( /* find NNEIGH neighbors for pixel */
319     short nl[NNEIGH][2],
320     int nd[NNEIGH],
321     int h,
322     int v,
323     register short (*rnl)[NNEIGH]
324     )
325 gwlarson 3.2 {
326     int nn = 0;
327 gwlarson 3.9 int d, n, hoff;
328 gwlarson 3.2 register int h2, n2;
329    
330 gwlarson 3.3 nd[NNEIGH-1] = MAXRAD2;
331 gwlarson 3.9 for (hoff = 0; hoff < hres; hoff = (hoff<=0) - hoff) {
332 gwlarson 3.2 h2 = h + hoff;
333 schorsch 3.14 if ((h2 < 0) | (h2 >= hres))
334 gwlarson 3.2 continue;
335 gwlarson 3.3 if ((h2-h)*(h2-h) >= nd[NNEIGH-1])
336 gwlarson 3.2 break;
337     for (n = 0; n < NNEIGH && rnl[h2][n] < NINF; n++) {
338     d = (h2-h)*(h2-h) + (v-rnl[h2][n])*(v-rnl[h2][n]);
339 schorsch 3.14 if ((d == 0) | (d >= nd[NNEIGH-1]))
340 gwlarson 3.2 continue;
341 gwlarson 3.3 if (nn < NNEIGH) /* insert neighbor */
342     nn++;
343     for (n2 = nn; n2--; ) {
344 gwlarson 3.2 if (!n2 || d >= nd[n2-1]) {
345     nd[n2] = d;
346     nl[n2][0] = h2;
347     nl[n2][1] = rnl[h2][n];
348     break;
349     }
350     nd[n2] = nd[n2-1];
351     nl[n2][0] = nl[n2-1][0];
352     nl[n2][1] = nl[n2-1][1];
353     }
354     }
355     }
356     return(nn);
357     }
358    
359    
360 schorsch 3.15 static void
361     meet_neighbors( /* run through samples and their neighbors */
362     int occ,
363     sampfunc *nf,
364     double *dp
365     )
366 gwlarson 3.2 {
367     short ln[NNEIGH][2];
368 gwlarson 3.9 int nd[NNEIGH];
369 gwlarson 3.2 int h, v, n, v2;
370     register short (*rnl)[NNEIGH];
371     /* initialize bottom row list */
372     rnl = (short (*)[NNEIGH])malloc(NNEIGH*sizeof(short)*hres);
373     CHECK(rnl==NULL, SYSTEM, "out of memory in meet_neighbors");
374     for (h = 0; h < hres; h++) {
375     for (n = v = 0; v < vres; v++)
376     if (CHK4(pixFlags, v*hres+h)) {
377     rnl[h][n++] = v;
378     if (n >= NNEIGH)
379     break;
380     }
381     while (n < NNEIGH)
382     rnl[h][n++] = NINF;
383     }
384     v = 0; /* do each row */
385     for ( ; ; ) {
386     for (h = 0; h < hres; h++) {
387 gwlarson 3.9 if (!CHK4(pixFlags, v*hres+h) != !occ)
388     continue; /* occupancy mismatch */
389     /* find neighbors */
390     n = findneigh(ln, nd, h, v, rnl);
391     /* call on neighbors */
392     (*nf)(h, v, ln, nd, n, dp);
393 gwlarson 3.2 }
394     if (++v >= vres) /* reinitialize row list */
395     break;
396     for (h = 0; h < hres; h++)
397     for (v2 = rnl[h][NNEIGH-1]+1; v2 < vres; v2++) {
398     if (v2 - v > v - rnl[h][0])
399     break; /* not close enough */
400     if (CHK4(pixFlags, v2*hres+h)) {
401     for (n = 0; n < NNEIGH-1; n++)
402     rnl[h][n] = rnl[h][n+1];
403     rnl[h][NNEIGH-1] = v2;
404     }
405     }
406     }
407 greg 3.10 free((void *)rnl); /* free row list */
408 gwlarson 3.1 }