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.8 by gwlarson, Tue Mar 9 15:08:22 1999 UTC

# Line 10 | Line 10 | static char SCCSid[] = "$SunId$ SGI";
10  
11   #include "holo.h"
12   #include "view.h"
13 + #include "random.h"
14  
15   #ifndef DEPS
16 < #define DEPS            0.01    /* depth epsilon */
16 > #define DEPS            0.02    /* depth epsilon */
17   #endif
18 + #ifndef PEPS
19 + #define PEPS            0.04    /* pixel value epsilon */
20 + #endif
21   #ifndef MAXRAD
22   #define MAXRAD          64      /* maximum kernel radius */
23   #endif
# Line 39 | Line 43 | static char SCCSid[] = "$SunId$ SGI";
43  
44   static int4     *pixFlags;      /* pixel occupancy flags */
45   static float    pixWeight[MAXRAD2];     /* pixel weighting function */
46 + static short    isqrttab[MAXRAD2];      /* integer square root table */
47  
48 + #define isqrt(i2)       ((int)isqrttab[(int)(i2)])
49 +
50   extern VIEW     myview;         /* current output view */
51   extern COLOR    *mypixel;       /* pixels being rendered */
52   extern float    *myweight;      /* weights (used to compute final pixels) */
# Line 107 | Line 114 | register HDBEAMI       *hb;
114   int
115   kill_occl(h, v, nl, n)          /* check for occlusion errors */
116   int     h, v;
117 < short   nl[NNEIGH][2];
117 > register short  nl[NNEIGH][2];
118   int     n;
119   {
120          short   forequad[2][2];
121          int     d;
122 <        register int4   i;
122 >        register int4   i, p;
123  
124 <        if (n <= 0)
124 >        if (n <= 0) {
125 > #ifdef DEBUG
126 >                error(WARNING, "neighborless sample in kill_occl");
127 > #endif
128                  return(1);
129 +        }
130 +        p = v*hres + h;
131          forequad[0][0] = forequad[0][1] = forequad[1][0] = forequad[1][1] = 0;
132          for (i = n; i--; ) {
133                  d = (h-nl[i][0])*(h-nl[i][0]) + (v-nl[i][1])*(v-nl[i][1]);
134 <                if (mydepth[nl[i][1]*hres+nl[i][0]] <
135 <                                mydepth[v*hres+h]*(1.-DEPS*d))
134 >                d = isqrt(d);
135 >                if (mydepth[nl[i][1]*hres+nl[i][0]]*(1.+DEPS*d) < mydepth[p])
136                          forequad[nl[i][0]<h][nl[i][1]<v] = 1;
137          }
138 <        if (forequad[0][0]+forequad[0][1]+forequad[1][0]+forequad[1][1] > 1) {
139 <                i = v*hres + h;
140 <                setcolor(mypixel[i], 0., 0., 0.);
129 <                myweight[i] = 0.;       /* occupancy reset afterwards */
138 >        if (forequad[0][0]+forequad[0][1]+forequad[1][0]+forequad[1][1] > 2) {
139 >                setcolor(mypixel[p], 0., 0., 0.);
140 >                myweight[p] = 0.;       /* occupancy reset afterwards */
141          }
142          return(1);
143   }
144  
145  
146   int
147 < grow_samp(h, v, nl, n)          /* grow sample point appropriately */
147 > smooth_samp(h, v, nl, n)        /* grow sample point smoothly */
148   int     h, v;
149   register short  nl[NNEIGH][2];
150   int     n;
151   {
152 +        int     dis[NNEIGH], ndis;
153          COLOR   mykern[MAXRAD2];
142 –        float   mykw[MAXRAD2];
154          int4    maxr2;
155 <        double  w;
155 >        double  d;
156          register int4   p, r2;
157 <        int     maxr, h2, v2;
157 >        int     i, r, maxr, h2, v2;
158  
159          if (n <= 0)
160                  return(1);
161          p = v*hres + h;                         /* build kernel values */
162          maxr2 = (h-nl[n-1][0])*(h-nl[n-1][0]) + (v-nl[n-1][1])*(v-nl[n-1][1]);
163          DCHECK(maxr2>=MAXRAD2, CONSISTENCY, "out of range neighbor");
164 <        for (r2 = maxr2+1; --r2; ) {
165 <                copycolor(mykern[r2], mypixel[p]);
166 <                mykw[r2] = pixWeight[r2];
167 <                if (2*r2 >= maxr2)              /* soften skirt */
168 <                        mykw[r2] *= (2*(maxr2-r2)+1.0)/maxr2;
169 <                scalecolor(mykern[r2], mykw[r2]);
164 >        maxr = isqrt(maxr2);
165 >        for (v2 = 1; v2 <= maxr; v2++)
166 >                for (h2 = 0; h2 <= v2; h2++) {
167 >                        r2 = h2*h2 + v2*v2;
168 >                        if (r2 > maxr2) break;
169 >                        copycolor(mykern[r2], mypixel[p]);
170 >                        scalecolor(mykern[r2], pixWeight[r2]);
171 >                }
172 >        ndis = 0;                               /* find discontinuities */
173 >        for (i = n; i--; ) {
174 >                r2 = (h-nl[i][0])*(h-nl[i][0]) + (v-nl[i][1])*(v-nl[i][1]);
175 >                r = isqrt(r2);
176 >                d = mydepth[nl[i][1]*hres+nl[i][0]] / mydepth[p];
177 >                d = d>=1. ? d-1. : 1.-d;
178 >                if (d > r*DEPS || bigdiff(mypixel[p],
179 >                                mypixel[nl[i][1]*hres+nl[i][0]], r*PEPS))
180 >                        dis[ndis++] = i;
181          }
182 <        maxr = sqrt((double)maxr2) + .99;       /* stamp out that kernel */
182 >                                                /* stamp out that kernel */
183          for (v2 = v-maxr; v2 <= v+maxr; v2++) {
184 <                if (v2 < 0) continue;
185 <                if (v2 >= vres) break;
184 >                if (v2 < 0) v2 = 0;
185 >                else if (v2 >= vres) break;
186                  for (h2 = h-maxr; h2 <= h+maxr; h2++) {
187 <                        if (h2 < 0) continue;
188 <                        if (h2 >= hres) break;
189 <                        r2 = (v2-v)*(v2-v) + (h2-h)*(h2-h);
187 >                        if (h2 < 0) h2 = 0;
188 >                        else if (h2 >= hres) break;
189 >                        r2 = (h2-h)*(h2-h) + (v2-v)*(v2-v);
190                          if (r2 > maxr2) continue;
191                          if (CHK4(pixFlags, v2*hres+h2))
192                                  continue;       /* occupied */
193 +                        for (i = ndis; i--; ) {
194 +                                r = (h2-nl[dis[i]][0])*(h2-nl[dis[i]][0]) +
195 +                                        (v2-nl[dis[i]][1])*(v2-nl[dis[i]][1]);
196 +                                if (r < r2) break;
197 +                        }
198 +                        if (i >= 0) continue;   /* outside edge */
199                          addcolor(mypixel[v2*hres+h2], mykern[r2]);
200 <                        myweight[v2*hres+h2] += mykw[r2]*myweight[v*hres+h];
200 >                        myweight[v2*hres+h2] += pixWeight[r2] * myweight[p];
201                  }
202          }
203          return(1);
204   }
205  
206  
207 < pixFlush()                      /* done with beams -- flush pixel values */
207 > int
208 > random_samp(h, v, nl, n, rf)    /* grow sample point noisily */
209 > int     h, v;
210 > register short  nl[NNEIGH][2];
211 > int     n;
212 > double  *rf;
213   {
214 +        double  tv = *rf * (1.-1./NNEIGH);
215 +        int4    maxr2;
216 +        register int4   p, r2;
217 +        register int    i;
218 +        int     maxr, h2, v2;
219 +        COLOR   ctmp;
220 +
221 +        if (n <= 0)
222 +                return(1);
223 +        p = v*hres + h;                         /* compute kernel radius */
224 +        maxr2 = (h-nl[n-1][0])*(h-nl[n-1][0]) + (v-nl[n-1][1])*(v-nl[n-1][1]);
225 +        DCHECK(maxr2>=MAXRAD2, CONSISTENCY, "out of range neighbor");
226 +        maxr = isqrt(maxr2);
227 +                                                /* sample kernel */
228 +        for (v2 = v-maxr; v2 <= v+maxr; v2++) {
229 +                if (v2 < 0) v2 = 0;
230 +                else if (v2 >= vres) break;
231 +                for (h2 = h-maxr; h2 <= h+maxr; h2++) {
232 +                        if (h2 < 0) h2 = 0;
233 +                        else if (h2 >= hres) break;
234 +                        r2 = (h2-h)*(h2-h) + (v2-v)*(v2-v);
235 +                        if (r2 > maxr2) continue;
236 +                        if (CHK4(pixFlags, v2*hres+h2))
237 +                                continue;       /* occupied */
238 +                        if (frandom() < tv) {   /* use neighbor instead? */
239 +                                i = random() % n;
240 +                                r2 = nl[i][1]*hres + nl[i][0];
241 +                                copycolor(ctmp, mypixel[r2]);
242 +                                r2 = (h2-nl[i][0])*(h2-nl[i][0]) +
243 +                                                (v2-nl[i][1])*(v2-nl[i][1]);
244 +                                if (r2 < MAXRAD2) {
245 +                                        scalecolor(ctmp, pixWeight[r2]);
246 +                                        addcolor(mypixel[v2*hres+h2], ctmp);
247 +                                        myweight[v2*hres+h2] += pixWeight[r2] *
248 +                                              myweight[nl[i][1]*hres+nl[i][0]];
249 +                                        continue;
250 +                                }
251 +                        }
252 +                                                /* use central sample */
253 +                        copycolor(ctmp, mypixel[p]);
254 +                        scalecolor(ctmp, pixWeight[r2]);
255 +                        addcolor(mypixel[v2*hres+h2], ctmp);
256 +                        myweight[v2*hres+h2] += pixWeight[r2] * myweight[p];
257 +                }
258 +        }
259 +        return(1);
260 + }
261 +
262 +
263 + pixFinish(ransamp)              /* done with beams -- compute pixel values */
264 + double  ransamp;
265 + {
266 +        if (pixWeight[0] <= FTINY)
267 +                init_wfunc();           /* initialize weighting function */
268          reset_flags();                  /* set occupancy flags */
269 <        meet_neighbors(kill_occl);      /* eliminate occlusion errors */
269 >        meet_neighbors(kill_occl,NULL); /* eliminate occlusion errors */
270          reset_flags();                  /* reset occupancy flags */
271 <        if (pixWeight[0] <= FTINY) {    /* initialize weighting function */
272 <                register int    r;
273 <                for (r = MAXRAD2; --r; )
274 <                        pixWeight[r] = G0NORM/sqrt((double)r);
188 <                pixWeight[0] = 1.;
189 <        }
190 <        meet_neighbors(grow_samp);      /* grow valid samples over image */
271 >        if (ransamp >= 0.)              /* spread samples over image */
272 >                meet_neighbors(random_samp,&ransamp);
273 >        else
274 >                meet_neighbors(smooth_samp,NULL);
275          free((char *)pixFlags);         /* free pixel flags */
276          pixFlags = NULL;
277   }
# Line 195 | Line 279 | pixFlush()                     /* done with beams -- flush pixel values
279  
280   reset_flags()                   /* allocate/set/reset occupancy flags */
281   {
282 <        register int    p;
282 >        register int4   p;
283  
284          if (pixFlags == NULL) {
285                  pixFlags = (int4 *)calloc(FL4NELS(hres*vres), sizeof(int4));
# Line 208 | Line 292 | reset_flags()                  /* allocate/set/reset occupancy flags
292   }
293  
294  
295 + init_wfunc()                    /* initialize weighting function */
296 + {
297 +        register int4   r2;
298 +        register double d;
299 +
300 +        for (r2 = MAXRAD2; --r2; ) {
301 +                d = sqrt((double)r2);
302 +                pixWeight[r2] = G0NORM/d;
303 +                isqrttab[r2] = d + 0.99;
304 +        }
305 +        pixWeight[0] = 1.;
306 +        isqrttab[0] = 0;
307 + }
308 +
309 +
310   int
311   findneigh(nl, h, v, rnl)        /* find NNEIGH neighbors for pixel */
312   short   nl[NNEIGH][2];
# Line 215 | Line 314 | int    h, v;
314   register short  (*rnl)[NNEIGH];
315   {
316          int     nn = 0;
317 <        int4    d, ld, nd[NNEIGH+1];
317 >        int4    d, nd[NNEIGH];
318          int     n, hoff;
319          register int    h2, n2;
320  
321 <        ld = MAXRAD2;
321 >        nd[NNEIGH-1] = MAXRAD2;
322          for (hoff = 1; hoff < hres; hoff = (hoff<0) - hoff) {
323                  h2 = h + hoff;
324                  if (h2 < 0 | h2 >= hres)
325                          continue;
326 <                if ((h2-h)*(h2-h) >= ld)
326 >                if ((h2-h)*(h2-h) >= nd[NNEIGH-1])
327                          break;
328                  for (n = 0; n < NNEIGH && rnl[h2][n] < NINF; n++) {
329                          d = (h2-h)*(h2-h) + (v-rnl[h2][n])*(v-rnl[h2][n]);
330 <                        if (d >= ld)
330 >                        if (d >= nd[NNEIGH-1])
331                                  continue;
332 <                        for (n2 = nn; ; n2--)   {       /* insert neighbor */
332 >                        if (nn < NNEIGH)        /* insert neighbor */
333 >                                nn++;
334 >                        for (n2 = nn; n2--; )   {
335                                  if (!n2 || d >= nd[n2-1]) {
336                                          nd[n2] = d;
337                                          nl[n2][0] = h2;
# Line 241 | Line 342 | register short (*rnl)[NNEIGH];
342                                  nl[n2][0] = nl[n2-1][0];
343                                  nl[n2][1] = nl[n2-1][1];
344                          }
244 –                        if (nn < NNEIGH)
245 –                                nn++;
246 –                        else
247 –                                ld = nd[NNEIGH-1];
345                  }
346          }
347          return(nn);
348   }
349  
350  
351 < meet_neighbors(nf)              /* run through samples and their neighbors */
351 > meet_neighbors(nf, dp)          /* run through samples and their neighbors */
352   int     (*nf)();
353 + char    *dp;
354   {
355          short   ln[NNEIGH][2];
356          int     h, v, n, v2;
# Line 276 | Line 374 | int    (*nf)();
374                          if (!CHK4(pixFlags, v*hres+h))
375                                  continue;       /* no one home */
376                          n = findneigh(ln, h, v, rnl);
377 <                        (*nf)(h, v, ln, n);     /* call on neighbors */
377 >                        (*nf)(h, v, ln, n, dp); /* call on neighbors */
378                  }
379                  if (++v >= vres)                /* reinitialize row list */
380                          break;

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines