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.1 by greg, Fri Oct 19 04:14:29 2012 UTC vs.
Revision 2.6 by greg, Wed Sep 25 05:03:10 2013 UTC

# Line 51 | Line 51 | add_bsdf_data(double theta_out, double phi_out, double
51          if (!isDSF)
52                  val *= ovec[2];         /* convert from BSDF to DSF */
53  
54 +                                        /* update BSDF histogram */
55 +        if (val < BSDF2BIG*ovec[2] && val > BSDF2SML*ovec[2])
56 +                ++bsdf_hist[histndx(val/ovec[2])];
57 +
58          pos_from_vec(pos, ovec);
59  
60          dsf_grid[pos[0]][pos[1]].vsum += val;
# Line 185 | Line 189 | cull_values(void)
189                  }
190   }
191  
192 + /* Compute minimum BSDF from histogram and clear it */
193 + static void
194 + comp_bsdf_min()
195 + {
196 +        int     cnt;
197 +        int     i, target;
198 +
199 +        cnt = 0;
200 +        for (i = HISTLEN; i--; )
201 +                cnt += bsdf_hist[i];
202 +        if (!cnt) {                             /* shouldn't happen */
203 +                bsdf_min = 0;
204 +                return;
205 +        }
206 +        target = cnt/100;                       /* ignore bottom 1% */
207 +        cnt = 0;
208 +        for (i = 0; cnt <= target; i++)
209 +                cnt += bsdf_hist[i];
210 +        bsdf_min = histval(i-1);
211 +        memset(bsdf_hist, 0, sizeof(bsdf_hist));
212 + }
213 +
214 + /* Find n nearest sub-sampled neighbors to the given grid position */
215 + static int
216 + get_neighbors(int neigh[][2], int n, const int i, const int j)
217 + {
218 +        int     k = 0;
219 +        int     r;
220 +                                                /* search concentric squares */
221 +        for (r = 1; r < GRIDRES; r++) {
222 +                int     ii, jj;
223 +                for (ii = i-r; ii <= i+r; ii++) {
224 +                        int     jstep = 1;
225 +                        if (ii < 0) continue;
226 +                        if (ii >= GRIDRES) break;
227 +                        if ((i-r < ii) & (ii < i+r))
228 +                                jstep = r<<1;
229 +                        for (jj = j-r; jj <= j+r; jj += jstep) {
230 +                                if (jj < 0) continue;
231 +                                if (jj >= GRIDRES) break;
232 +                                if (dsf_grid[ii][jj].nval) {
233 +                                        neigh[k][0] = ii;
234 +                                        neigh[k][1] = jj;
235 +                                        if (++k >= n)
236 +                                                return(n);
237 +                                }
238 +                        }
239 +                }
240 +        }
241 +        return(k);
242 + }
243 +
244 + /* Adjust coded radius for the given grid position based on neighborhood */
245 + static int
246 + adj_coded_radius(const int i, const int j)
247 + {
248 +        const double    max_frac = 0.33;
249 +        const double    rad0 = R2ANG(dsf_grid[i][j].crad);
250 +        double          currad = RSCA * rad0;
251 +        int             neigh[5][2];
252 +        int             n;
253 +        FVECT           our_dir;
254 +
255 +        ovec_from_pos(our_dir, i, j);
256 +        n = get_neighbors(neigh, 5, i, j);
257 +        while (n--) {
258 +                FVECT   their_dir;
259 +                double  max_ratio, rad_ok2;
260 +                                                /* check our value at neighbor */
261 +                ovec_from_pos(their_dir, neigh[n][0], neigh[n][1]);
262 +                max_ratio = max_frac * dsf_grid[neigh[n][0]][neigh[n][1]].vsum
263 +                                / dsf_grid[i][j].vsum;
264 +                if (max_ratio >= 1)
265 +                        continue;
266 +                rad_ok2 = (DOT(their_dir,our_dir) - 1.)/log(max_ratio);
267 +                if (rad_ok2 >= currad*currad)
268 +                        continue;               /* value fraction OK */
269 +                currad = sqrt(rad_ok2);         /* else reduce lobe radius */
270 +                if (currad <= rad0)             /* limit how small we'll go */
271 +                        return(dsf_grid[i][j].crad);
272 +        }
273 +        return(ANG2R(currad));                  /* encode selected radius */
274 + }
275 +
276   /* Count up filled nodes and build RBF representation from current grid */
277   RBFNODE *
278   make_rbfrep(void)
# Line 193 | Line 281 | make_rbfrep(void)
281          double  lastVar, thisVar = 100.;
282          int     nn;
283          RBFNODE *newnode;
284 +        RBFVAL  *itera;
285          int     i, j;
286                                  /* compute RBF radii */
287          compute_radii();
# Line 202 | Line 291 | make_rbfrep(void)
291          for (i = 0; i < GRIDRES; i++)
292              for (j = 0; j < GRIDRES; j++)
293                  nn += dsf_grid[i][j].nval;
294 +                                /* compute minimum BSDF */
295 +        comp_bsdf_min();
296                                  /* allocate RBF array */
297          newnode = (RBFNODE *)malloc(sizeof(RBFNODE) + sizeof(RBFVAL)*(nn-1));
298 <        if (newnode == NULL) {
299 <                fprintf(stderr, "%s: Out of memory in make_rbfrep()\n", progname);
209 <                exit(1);
210 <        }
298 >        if (newnode == NULL)
299 >                goto memerr;
300          newnode->ord = -1;
301          newnode->next = NULL;
302          newnode->ejl = NULL;
# Line 222 | Line 311 | make_rbfrep(void)
311              for (j = 0; j < GRIDRES; j++)
312                  if (dsf_grid[i][j].nval) {
313                          newnode->rbfa[nn].peak = dsf_grid[i][j].vsum;
314 <                        newnode->rbfa[nn].crad = RSCA*dsf_grid[i][j].crad + .5;
314 >                        newnode->rbfa[nn].crad = adj_coded_radius(i, j);
315                          newnode->rbfa[nn].gx = i;
316                          newnode->rbfa[nn].gy = j;
317                          ++nn;
318                  }
319                                  /* iterate to improve interpolation accuracy */
320 +        itera = (RBFVAL *)malloc(sizeof(RBFVAL)*newnode->nrbf);
321 +        if (itera == NULL)
322 +                goto memerr;
323 +        memcpy(itera, newnode->rbfa, sizeof(RBFVAL)*newnode->nrbf);
324          do {
325                  double  dsum = 0, dsum2 = 0;
326                  nn = 0;
# Line 237 | Line 330 | make_rbfrep(void)
330                                  FVECT   odir;
331                                  double  corr;
332                                  ovec_from_pos(odir, i, j);
333 <                                newnode->rbfa[nn++].peak *= corr =
333 >                                itera[nn++].peak *= corr =
334                                          dsf_grid[i][j].vsum /
335                                                  eval_rbfrep(newnode, odir);
336 <                                dsum += corr - 1.;
337 <                                dsum2 += (corr-1.)*(corr-1.);
336 >                                dsum += 1. - corr;
337 >                                dsum2 += (1.-corr)*(1.-corr);
338                          }
339 +                memcpy(newnode->rbfa, itera, sizeof(RBFVAL)*newnode->nrbf);
340                  lastVar = thisVar;
341                  thisVar = dsum2/(double)nn;
342   #ifdef DEBUG
# Line 252 | Line 346 | make_rbfrep(void)
346   #endif
347          } while (--niter > 0 && lastVar-thisVar > 0.02*lastVar);
348  
349 +        free(itera);
350          nn = 0;                 /* compute sum for normalization */
351          while (nn < newnode->nrbf)
352                  newnode->vtotal += rbf_volume(&newnode->rbfa[nn++]);
353 <
353 > #ifdef DEBUG
354 >        fprintf(stderr, "Integrated DSF at (%.1f,%.1f) deg. is %.2f\n",
355 >                        get_theta180(newnode->invec), get_phi360(newnode->invec),
356 >                        newnode->vtotal);
357 > #endif
358          insert_dsf(newnode);
359  
360          return(newnode);
361 + memerr:
362 +        fprintf(stderr, "%s: Out of memory in make_rbfrep()\n", progname);
363 +        exit(1);
364   }

Diff Legend

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