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

Comparing ray/src/cv/bsdfrbf.c (file contents):
Revision 2.10 by greg, Fri Oct 18 02:49:30 2013 UTC vs.
Revision 2.11 by greg, Sat Oct 19 00:11:50 2013 UTC

# Line 15 | Line 15 | static const char RCSid[] = "$Id$";
15   #include "bsdfrep.h"
16  
17   #ifndef MINRSCA
18 < #define MINRSCA         0.5             /* minimum radius scaling factor */
18 > #define MINRSCA         1.0             /* minimum radius scaling factor */
19   #endif
20   #ifndef MAXRSCA
21   #define MAXRSCA         2.7             /* maximum radius scaling factor */
22   #endif
23 < #ifndef DIFFTHRESH
24 < #define DIFFTHRESH      0.2             /* culling difference threshold */
23 > #ifndef VARTHRESH
24 > #define VARTHRESH       0.0015          /* culling variance threshold */
25   #endif
26 + #ifndef DIFFMAX2
27 + #define DIFFMAX2        (16.*VARTHRESH) /* maximum ignored sample variance */
28 + #endif
29   #ifndef MAXFRAC
30   #define MAXFRAC         0.5             /* maximum contribution to neighbor */
31   #endif
# Line 75 | Line 78 | add_bsdf_data(double theta_out, double phi_out, double
78          dsf_grid[pos[0]][pos[1]].nval++;
79   }
80  
78 /* Check if the two DSF values are significantly different */
79 static int
80 big_diff(double ref, double tst)
81 {
82        if (ref > 0) {
83                tst = tst/ref - 1.;
84                if (tst < 0) tst = -tst;
85        } else
86                tst *= 50.;
87        return(tst > DIFFTHRESH);
88 }
89
81   /* Compute radii for non-empty bins */
82   /* (distance to furthest empty bin for which non-empty test bin is closest) */
83   static void
# Line 139 | Line 130 | compute_radii(void)
130              }
131          for (i = 0; i < GRIDRES; i++)           /* grow radii where uniform */
132              for (j = 0; j < GRIDRES; j++) {
133 <                double  midmean;
133 >                double  midmean = 0.0;
134 >                int     nsum = 0;
135                  if (!dsf_grid[i][j].nval)
136                          continue;
137 <                midmean = dsf_grid[i][j].vsum / (double)dsf_grid[i][j].nval;
138 <                r = R2ANG(dsf_grid[i][j].crad)*(MAXRSCA*GRIDRES/M_PI);
137 >                r = 1;                          /* avg. immediate neighbors */
138 >                for (ii = i-r; ii <= i+r; ii++) {
139 >                    if (ii < 0) continue;
140 >                    if (ii >= GRIDRES) break;
141 >                    for (jj = j-r; jj <= j+r; jj++) {
142 >                        if (jj < 0) continue;
143 >                        if (jj >= GRIDRES) break;
144 >                        midmean += dsf_grid[ii][jj].vsum;
145 >                        nsum += dsf_grid[ii][jj].nval;
146 >                    }
147 >                }
148 >                midmean /= (double)nsum;
149                  while (++r < GRIDRES) {         /* attempt to grow perimeter */
150 +                    double      diff2sum = 0.0;
151 +                    nsum = 0;
152                      for (ii = i-r; ii <= i+r; ii++) {
153                          int     jstep = 1;
154                          if (ii < 0) continue;
# Line 152 | Line 156 | compute_radii(void)
156                          if ((i-r < ii) & (ii < i+r))
157                                  jstep = r<<1;
158                          for (jj = j-r; jj <= j+r; jj += jstep) {
159 +                                double  d2;
160                                  if (jj < 0) continue;
161                                  if (jj >= GRIDRES) break;
162 <                                if (dsf_grid[ii][jj].nval && big_diff(midmean,
163 <                                                dsf_grid[ii][jj].vsum /
164 <                                                (double)dsf_grid[ii][jj].nval))
165 <                                        goto hit_diff;
162 >                                if (!dsf_grid[ii][jj].nval)
163 >                                        continue;
164 >                                d2 = midmean - dsf_grid[ii][jj].vsum /
165 >                                                (double)dsf_grid[ii][jj].nval;
166 >                                d2 *= d2;
167 >                                if (d2 > DIFFMAX2*midmean*midmean)
168 >                                        goto escape;
169 >                                diff2sum += d2;
170 >                                ++nsum;
171                          }
172                      }
173 +                    if (diff2sum > VARTHRESH*midmean*midmean*(double)nsum)
174 +                        break;
175                  }
176 < hit_diff:       --r;
177 <                dsf_grid[i][j].crad = ANG2R(r*(M_PI/MAXRSCA/GRIDRES));
178 <                if (dsf_grid[i][j].crad < cradmin)
179 <                        dsf_grid[i][j].crad = cradmin;
176 > escape:         --r;
177 >                r = ANG2R(r*(M_PI/MAXRSCA/GRIDRES));
178 >                if (r < cradmin)
179 >                        r = cradmin;
180 >                if (dsf_grid[i][j].crad < r)
181 >                        dsf_grid[i][j].crad = r;
182              }
183                                                  /* blur radii over hemisphere */
184          memset(fill_grid, 0, sizeof(fill_grid));
# Line 194 | Line 208 | hit_diff:      --r;
208                          dsf_grid[i][j].crad = fill_grid[i][j]/fill_cnt[i][j];
209   }
210  
211 + /* Radius comparison for qsort() */
212 + static int
213 + radius_cmp(const void *p1, const void *p2)
214 + {
215 +        return( (int)dsf_grid[0][*(const int *)p1].crad -
216 +                (int)dsf_grid[0][*(const int *)p2].crad );
217 + }
218 +
219   /* Cull points for more uniform distribution, leave all nval 0 or 1 */
220   static void
221   cull_values(void)
222   {
223 +        int     indx[GRIDRES*GRIDRES];
224          FVECT   ovec0, ovec1;
225          double  maxang, maxang2;
226 <        int     i, j, ii, jj, r;
226 >        int     i, j, k, ii, jj, r;
227 >                                                /* sort by radius first */
228 >        for (k = GRIDRES*GRIDRES; k--; )
229 >                indx[k] = k;
230 >        qsort(indx, GRIDRES*GRIDRES, sizeof(int), &radius_cmp);
231                                                  /* simple greedy algorithm */
232 <        for (i = 0; i < GRIDRES; i++)
233 <            for (j = 0; j < GRIDRES; j++) {
232 >        for (k = GRIDRES*GRIDRES; k--; ) {
233 >                i = indx[k]/GRIDRES;            /* from biggest radius down */
234 >                j = indx[k] - i*GRIDRES;
235                  if (!dsf_grid[i][j].nval)
236                          continue;
237                  if (!dsf_grid[i][j].crad)
238 <                        continue;               /* shouldn't happen */
238 >                        break;          /* shouldn't happen */
239                  ovec_from_pos(ovec0, i, j);
240                  maxang = 2.*R2ANG(dsf_grid[i][j].crad);
241                                                  /* clamp near horizon */
# Line 236 | Line 264 | cull_values(void)
264                          dsf_grid[ii][jj].nval = 0;
265                      }
266                  }
267 <            }
267 >        }
268                                                  /* final averaging pass */
269          for (i = 0; i < GRIDRES; i++)
270              for (j = 0; j < GRIDRES; j++)
# Line 246 | Line 274 | cull_values(void)
274                  }
275   }
276  
277 < /* Compute minimum BSDF from histogram and clear it */
277 > /* Compute minimum BSDF from histogram (does not clear) */
278   static void
279   comp_bsdf_min()
280   {
# Line 265 | Line 293 | comp_bsdf_min()
293          for (i = 0; cnt <= target; i++)
294                  cnt += bsdf_hist[i];
295          bsdf_min = histval(i-1);
268        memset(bsdf_hist, 0, sizeof(bsdf_hist));
296   }
297  
298   /* Find n nearest sub-sampled neighbors to the given grid position */

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines