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

Comparing ray/src/hd/rhpict2.c (file contents):
Revision 3.2 by gwlarson, Fri Mar 5 17:03:38 1999 UTC vs.
Revision 3.15 by schorsch, Thu Jan 1 11:21:55 2004 UTC

# Line 1 | Line 1
1 /* Copyright (c) 1999 Silicon Graphics, Inc. */
2
1   #ifndef lint
2 < static char SCCSid[] = "$SunId$ SGI";
2 > static const char       RCSid[] = "$Id$";
3   #endif
6
4   /*
5   * Rendering routines for rhpict.
6   */
7  
8 + #include <string.h>
9 +
10   #include "holo.h"
11   #include "view.h"
12 + #include "random.h"
13  
14   #ifndef DEPS
15 < #define DEPS            0.01    /* depth epsilon */
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 */
24 > #define NNEIGH          5       /* find this many neighbors */
25   #endif
26  
27   #define NINF            16382
# Line 34 | Line 37 | static char SCCSid[] = "$SunId$ SGI";
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))
40 > #define CLR4ALL(f,n)    memset((char *)(f),'\0',FL4NELS(n)*sizeof(int32))
41   #endif
42  
43 < static int4     *pixFlags;      /* pixel occupancy flags */
43 > static int32    *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)       (isqrttab[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 + extern void pixFinish(double ransamp);
56 + extern void pixBeam(BEAM *bp, HDBEAMI *hb);
57  
58 < pixBeam(bp, hb)                 /* render a particular beam */
59 < BEAM    *bp;
60 < register HDBEAMI        *hb;
58 > 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   {
76          GCOORD  gc[2];
77          register RAYVAL *rv;
# Line 57 | Line 79 | register HDBEAMI       *hb;
79          double  d, prox;
80          COLOR   col;
81          int     n;
82 <        register int4   p;
82 >        register int32  p;
83  
84          if (!hdbcoord(gc, hb->h, hb->b))
85                  error(CONSISTENCY, "bad beam in render_beam");
# Line 104 | Line 126 | register HDBEAMI       *hb;
126   }
127  
128  
129 < int
130 < kill_occl(h, v, nl, n)          /* check for occlusion errors */
131 < int     h, v;
132 < short   nl[NNEIGH][2];
133 < int     n;
129 > 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   {
139          short   forequad[2][2];
140          int     d;
141 <        register int4   i;
141 >        register int    i;
142 >        register int32  p;
143  
144 <        if (n <= 0)
144 >        if (n <= 0) {
145 > #ifdef DEBUG
146 >                error(WARNING, "neighborless sample in kill_occl");
147 > #endif
148                  return(1);
149 +        }
150 +        p = v*hres + h;
151          forequad[0][0] = forequad[0][1] = forequad[1][0] = forequad[1][1] = 0;
152          for (i = n; i--; ) {
153 <                d = (h-nl[i][0])*(h-nl[i][0]) + (v-nl[i][1])*(v-nl[i][1]);
154 <                if (mydepth[nl[i][1]*hres+nl[i][0]] <
123 <                                mydepth[v*hres+h]*(1.-DEPS*d))
153 >                d = isqrt(nd[i]);
154 >                if (mydepth[nl[i][1]*hres+nl[i][0]]*(1.+DEPS*d) < mydepth[p])
155                          forequad[nl[i][0]<h][nl[i][1]<v] = 1;
156          }
157 <        if (forequad[0][0]+forequad[0][1]+forequad[1][0]+forequad[1][1] > 1) {
158 <                i = v*hres + h;
159 <                setcolor(mypixel[i], 0., 0., 0.);
129 <                myweight[i] = 0.;       /* occupancy reset afterwards */
157 >        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          }
161          return(1);
162   }
163  
164  
165 < int
166 < grow_samp(h, v, nl, n)          /* grow sample point appropriately */
167 < int     h, v;
168 < register short  nl[NNEIGH][2];
169 < int     n;
165 > 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   {
175 +        int     dis[NNEIGH], ndis;
176          COLOR   mykern[MAXRAD2];
177 <        float   mykw[MAXRAD2];
178 <        int4    maxr2;
179 <        double  w;
180 <        register int4   p, r2;
181 <        int     maxr, h2, v2;
177 >        int     maxr2;
178 >        double  d;
179 >        register int32  p;
180 >        register int    r2;
181 >        int     i, r, maxr, h2, v2;
182  
183          if (n <= 0)
184                  return(1);
185          p = v*hres + h;                         /* build kernel values */
186 <        maxr2 = (h-nl[n-1][0])*(h-nl[n-1][0]) + (v-nl[n-1][1])*(v-nl[n-1][1]);
186 >        maxr2 = nd[n-1];
187          DCHECK(maxr2>=MAXRAD2, CONSISTENCY, "out of range neighbor");
188 <        for (r2 = maxr2+1; --r2; ) {
189 <                copycolor(mykern[r2], mypixel[p]);
190 <                mykw[r2] = pixWeight[r2];
191 <                if (2*r2 >= maxr2)              /* soften skirt */
192 <                        mykw[r2] *= (2*(maxr2-r2)+1.0)/maxr2;
193 <                scalecolor(mykern[r2], mykw[r2]);
188 >        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 >                r = isqrt(nd[i]);
199 >                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          }
205 <        maxr = sqrt((double)maxr2) + .99;       /* stamp out that kernel */
205 >                                                /* stamp out that kernel */
206          for (v2 = v-maxr; v2 <= v+maxr; v2++) {
207 <                if (v2 < 0) continue;
208 <                if (v2 >= vres) break;
207 >                if (v2 < 0) v2 = 0;
208 >                else if (v2 >= vres) break;
209                  for (h2 = h-maxr; h2 <= h+maxr; h2++) {
210 <                        if (h2 < 0) continue;
211 <                        if (h2 >= hres) break;
212 <                        r2 = (v2-v)*(v2-v) + (h2-h)*(h2-h);
210 >                        if (h2 < 0) h2 = 0;
211 >                        else if (h2 >= hres) break;
212 >                        r2 = (h2-h)*(h2-h) + (v2-v)*(v2-v);
213                          if (r2 > maxr2) continue;
214                          if (CHK4(pixFlags, v2*hres+h2))
215                                  continue;       /* occupied */
216 +                        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                          addcolor(mypixel[v2*hres+h2], mykern[r2]);
223 <                        myweight[v2*hres+h2] += mykw[r2]*myweight[v*hres+h];
223 >                        myweight[v2*hres+h2] += pixWeight[r2] * myweight[p];
224                  }
225          }
226          return(1);
227   }
228  
229  
230 < pixFlush()                      /* done with beams -- flush pixel values */
230 > 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   {
240 +        float   rnt[NNEIGH];
241 +        double  rvar;
242 +        register int32  p, pn;
243 +        register int    ni;
244 +
245 +        if (n <= 0)
246 +                return(1);
247 +        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 +        }
259 +        pn = nl[ni][1]*hres + nl[ni][0];
260 +        addcolor(mypixel[p], mypixel[pn]);
261 +        myweight[p] += myweight[pn];
262 +        return(1);
263 + }
264 +
265 +
266 + extern void
267 + pixFinish(              /* done with beams -- compute pixel values */
268 +        double  ransamp
269 + )
270 + {
271 +        if (pixWeight[0] <= FTINY)
272 +                init_wfunc();           /* initialize weighting function */
273          reset_flags();                  /* set occupancy flags */
274 <        meet_neighbors(kill_occl);      /* eliminate occlusion errors */
274 >        meet_neighbors(1,kill_occl,NULL); /* identify occlusion errors */
275          reset_flags();                  /* reset occupancy flags */
276 <        if (pixWeight[0] <= FTINY) {    /* initialize weighting function */
277 <                register int    r;
278 <                for (r = MAXRAD2; --r; )
279 <                        pixWeight[r] = G0NORM/sqrt((double)r);
280 <                pixWeight[0] = 1.;
189 <        }
190 <        meet_neighbors(grow_samp);      /* grow valid samples over image */
191 <        free((char *)pixFlags);         /* free pixel flags */
276 >        if (ransamp >= 0.)              /* spread samples over image */
277 >                meet_neighbors(0,random_samp,&ransamp);
278 >        else
279 >                meet_neighbors(1,smooth_samp,NULL);
280 >        free((void *)pixFlags);         /* free pixel flags */
281          pixFlags = NULL;
282   }
283  
284  
285 < reset_flags()                   /* allocate/set/reset occupancy flags */
285 > static void
286 > reset_flags(void)                       /* allocate/set/reset occupancy flags */
287   {
288 <        register int    p;
288 >        register int32  p;
289  
290          if (pixFlags == NULL) {
291 <                pixFlags = (int4 *)calloc(FL4NELS(hres*vres), sizeof(int4));
291 >                pixFlags = (int32 *)calloc(FL4NELS(hres*vres), sizeof(int32));
292                  CHECK(pixFlags==NULL, SYSTEM, "out of memory in reset_flags");
293          } else
294                  CLR4ALL(pixFlags, hres*vres);
# Line 208 | Line 298 | reset_flags()                  /* allocate/set/reset occupancy flags
298   }
299  
300  
301 < int
302 < findneigh(nl, h, v, rnl)        /* find NNEIGH neighbors for pixel */
213 < short   nl[NNEIGH][2];
214 < int     h, v;
215 < register short  (*rnl)[NNEIGH];
301 > static void
302 > init_wfunc(void)                        /* initialize weighting function */
303   {
304 +        register int    r2;
305 +        register double d;
306 +
307 +        for (r2 = MAXRAD2; --r2; ) {
308 +                d = sqrt((double)r2);
309 +                pixWeight[r2] = G0NORM/d;
310 +                isqrttab[r2] = d + 0.99;
311 +        }
312 +        pixWeight[0] = 1.;
313 +        isqrttab[0] = 0;
314 + }
315 +
316 +
317 + 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 + {
326          int     nn = 0;
327 <        int4    d, ld, nd[NNEIGH+1];
219 <        int     n, hoff;
327 >        int     d, n, hoff;
328          register int    h2, n2;
329  
330 <        ld = MAXRAD2;
331 <        for (hoff = 1; hoff < hres; hoff = (hoff<0) - hoff) {
330 >        nd[NNEIGH-1] = MAXRAD2;
331 >        for (hoff = 0; hoff < hres; hoff = (hoff<=0) - hoff) {
332                  h2 = h + hoff;
333 <                if (h2 < 0 | h2 >= hres)
333 >                if ((h2 < 0) | (h2 >= hres))
334                          continue;
335 <                if ((h2-h)*(h2-h) >= ld)
335 >                if ((h2-h)*(h2-h) >= nd[NNEIGH-1])
336                          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 <                        if (d >= ld)
339 >                        if ((d == 0) | (d >= nd[NNEIGH-1]))
340                                  continue;
341 <                        for (n2 = nn; ; n2--)   {       /* insert neighbor */
341 >                        if (nn < NNEIGH)        /* insert neighbor */
342 >                                nn++;
343 >                        for (n2 = nn; n2--; )   {
344                                  if (!n2 || d >= nd[n2-1]) {
345                                          nd[n2] = d;
346                                          nl[n2][0] = h2;
# Line 241 | Line 351 | register short (*rnl)[NNEIGH];
351                                  nl[n2][0] = nl[n2-1][0];
352                                  nl[n2][1] = nl[n2-1][1];
353                          }
244                        if (nn < NNEIGH)
245                                nn++;
246                        else
247                                ld = nd[NNEIGH-1];
354                  }
355          }
356          return(nn);
357   }
358  
359  
360 < meet_neighbors(nf)              /* run through samples and their neighbors */
361 < int     (*nf)();
360 > static void
361 > meet_neighbors( /* run through samples and their neighbors */
362 >        int     occ,
363 >        sampfunc *nf,
364 >        double  *dp
365 > )
366   {
367          short   ln[NNEIGH][2];
368 +        int     nd[NNEIGH];
369          int     h, v, n, v2;
370          register short  (*rnl)[NNEIGH];
371                                          /* initialize bottom row list */
# Line 273 | Line 384 | int    (*nf)();
384          v = 0;                          /* do each row */
385          for ( ; ; ) {
386                  for (h = 0; h < hres; h++) {
387 <                        if (!CHK4(pixFlags, v*hres+h))
388 <                                continue;       /* no one home */
389 <                        n = findneigh(ln, h, v, rnl);
390 <                        (*nf)(h, v, ln, n);     /* call on neighbors */
387 >                        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                  }
394                  if (++v >= vres)                /* reinitialize row list */
395                          break;
# Line 291 | Line 404 | int    (*nf)();
404                                  }
405                          }
406          }
407 <        free((char *)rnl);              /* free row list */
407 >        free((void *)rnl);              /* free row list */
408   }

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines