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

Comparing ray/src/rt/ambient.c (file contents):
Revision 2.74 by greg, Sat Apr 19 02:39:44 2014 UTC vs.
Revision 2.91 by greg, Thu Jun 19 16:26:55 2014 UTC

# Line 52 | Line 52 | static int  nunflshed = 0;     /* number of unflushed ambi
52   #endif
53  
54  
55 – static double  qambacc = 0.;            /* ambient accuracy to the 1/4 power */
55   static double  avsum = 0.;              /* computed ambient value sum (log) */
56   static unsigned int  navsum = 0;        /* number of values in avsum */
57   static unsigned int  nambvals = 0;      /* total number of indirect values */
# Line 110 | Line 109 | setambres(                             /* set ambient resolution */
109                                                  /* set min & max radii */
110          if (ar <= 0) {
111                  minarad = 0;
112 <                maxarad = thescene.cusize / 2.0;
112 >                maxarad = thescene.cusize*0.2;
113          } else {
114                  minarad = thescene.cusize / ar;
115 <                maxarad = 64 * minarad;                 /* heuristic */
116 <                if (maxarad > thescene.cusize / 2.0)
117 <                        maxarad = thescene.cusize / 2.0;
115 >                maxarad = 64.0 * minarad;               /* heuristic */
116 >                if (maxarad > thescene.cusize*0.2)
117 >                        maxarad = thescene.cusize*0.2;
118          }
119          if (minarad <= FTINY)
120 <                minarad = 10*FTINY;
120 >                minarad = 10.0*FTINY;
121          if (maxarad <= minarad)
122 <                maxarad = 64 * minarad;
122 >                maxarad = 64.0 * minarad;
123   }
124  
125  
# Line 129 | Line 128 | setambacc(                             /* set ambient accuracy */
128          double  newa
129   )
130   {
131 <        double  ambdiff;
132 <
133 <        if (newa < 0.0)
134 <                newa = 0.0;
135 <        ambdiff = fabs(newa - ambacc);
137 <        if (ambdiff >= .01 && (ambacc = newa) > FTINY) {
138 <                qambacc = sqrt(sqrt(ambacc));
131 >        static double   olda;           /* remember previous setting here */
132 >        
133 >        newa *= (newa > 0);
134 >        if (fabs(newa - olda) >= .05*(newa + olda)) {
135 >                ambacc = newa;
136                  if (nambvals > 0)
137                          sortambvals(1);         /* rebuild tree */
138          }
# Line 168 | Line 165 | setambient(void)                               /* initialize calculation */
165                  initambfile(0);                 /* file exists */
166                  lastpos = ftell(ambfp);
167                  while (readambval(&amb, ambfp))
168 <                        avinsert(avstore(&amb));
168 >                        avstore(&amb);
169                  nambshare = nambvals;           /* share loaded values */
170                  if (readonly) {
171                          sprintf(errmsg,
# Line 266 | Line 263 | ambnotify(                     /* record new modifier */
263  
264   /************ THE FOLLOWING ROUTINES DIFFER BETWEEN NEW & OLD ***************/
265  
266 < #ifdef NEWAMB
266 > #ifndef OLDAMB
267  
268   #define tfunc(lwr, x, upr)      (((x)-(lwr))/((upr)-(lwr)))
269  
270 + static int      plugaleak(RAY *r, AMBVAL *ap, FVECT anorm, double ang);
271   static double   sumambient(COLOR acol, RAY *r, FVECT rn, int al,
272                                  AMBTREE *at, FVECT c0, double s);
273   static int      makeambient(COLOR acol, RAY *r, FVECT rn, int al);
# Line 301 | Line 299 | multambient(           /* compute ambient component & multiply
299          if (ambacc <= FTINY) {                  /* no ambient storage */
300                  copycolor(acol, aval);
301                  rdepth++;
302 <                ok = doambient(acol, r, r->rweight, NULL, NULL, NULL, NULL);
302 >                ok = doambient(acol, r, r->rweight,
303 >                                NULL, NULL, NULL, NULL, NULL);
304                  rdepth--;
305                  if (!ok)
306                          goto dumbamb;
# Line 347 | Line 346 | dumbamb:                                       /* return global value */
346   }
347  
348  
349 < double
349 > /* Plug a potential leak where ambient cache value is occluded */
350 > static int
351 > plugaleak(RAY *r, AMBVAL *ap, FVECT anorm, double ang)
352 > {
353 >        const double    cost70sq = 0.1169778;   /* cos(70deg)^2 */
354 >        RAY             rtst;
355 >        FVECT           vdif;
356 >        double          normdot, ndotd, nadotd;
357 >        double          a, b, c, t[2];
358 >
359 >        ang += 2.*PI*(ang < 0);                 /* check direction flags */
360 >        if ( !(ap->corral>>(int)(ang*(16./PI)) & 1) )
361 >                return(0);
362 >        /*
363 >         * Generate test ray, targeting 20 degrees above sample point plane
364 >         * along surface normal from cache position.  This should be high
365 >         * enough to miss local geometry we don't really care about.
366 >         */
367 >        VSUB(vdif, ap->pos, r->rop);
368 >        normdot = DOT(anorm, r->ron);
369 >        ndotd = DOT(vdif, r->ron);
370 >        nadotd = DOT(vdif, anorm);
371 >        a = normdot*normdot - cost70sq;
372 >        b = 2.0*(normdot*ndotd - nadotd*cost70sq);
373 >        c = ndotd*ndotd - DOT(vdif,vdif)*cost70sq;
374 >        if (quadratic(t, a, b, c) != 2)
375 >                return(1);                      /* should rarely happen */
376 >        if (t[1] <= FTINY)
377 >                return(0);                      /* should fail behind test */
378 >        rayorigin(&rtst, SHADOW, r, NULL);
379 >        VSUM(rtst.rdir, vdif, anorm, t[1]);     /* further dist. > plane */
380 >        rtst.rmax = normalize(rtst.rdir);       /* short ray test */
381 >        while (localhit(&rtst, &thescene)) {    /* check for occluder */
382 >                if (rtst.ro->omod != OVOID &&
383 >                                (rtst.clipset == NULL ||
384 >                                        !inset(rtst.clipset, rtst.ro->omod)))
385 >                        return(1);              /* plug light leak */
386 >                VCOPY(rtst.rorg, rtst.rop);     /* skip invisible surface */
387 >                rtst.rmax -= rtst.rot;
388 >                rayclear(&rtst);
389 >        }
390 >        return(0);                              /* seems we're OK */
391 > }
392 >
393 >
394 > static double
395   sumambient(             /* get interpolated ambient value */
396          COLOR  acol,
397          RAY  *r,
# Line 357 | Line 401 | sumambient(            /* get interpolated ambient value */
401          FVECT  c0,
402          double  s
403   )
404 < {                                       /* initial limit is ambacc radians */
405 <        const double    maxangle = (ambacc-PI/2.)*pow(r->rweight,0.13) + PI/2.;
404 > {                       /* initial limit is 10 degrees plus ambacc radians */
405 >        const double    minangle = 10.0 * PI/180.;
406 >        double          maxangle = minangle + ambacc;
407          double          wsum = 0.0;
408          FVECT           ck0;
409          int             i, j;
410          AMBVAL          *av;
411 +
412 +        if (at->kid != NULL) {          /* sum children first */                                
413 +                s *= 0.5;
414 +                for (i = 0; i < 8; i++) {
415 +                        for (j = 0; j < 3; j++) {
416 +                                ck0[j] = c0[j];
417 +                                if (1<<j & i)
418 +                                        ck0[j] += s;
419 +                                if (r->rop[j] < ck0[j] - OCTSCALE*s)
420 +                                        break;
421 +                                if (r->rop[j] > ck0[j] + (1.0+OCTSCALE)*s)
422 +                                        break;
423 +                        }
424 +                        if (j == 3)
425 +                                wsum += sumambient(acol, r, rn, al,
426 +                                                        at->kid+i, ck0, s);
427 +                }
428 +                                        /* good enough? */
429 +                if (wsum >= 0.05 && s > minarad*10.0)
430 +                        return(wsum);
431 +        }
432 +                                        /* adjust maximum angle */
433 +        if (at->alist != NULL && (at->alist->lvl <= al) & (r->rweight < 0.6))
434 +                maxangle = (maxangle - PI/2.)*pow(r->rweight,0.13) + PI/2.;
435                                          /* sum this node */
436          for (av = at->alist; av != NULL; av = av->next) {
437 <                double  d, delta_r2, delta_t2;
437 >                double  u, v, d, delta_r2, delta_t2;
438                  COLOR   ct;
439                  FVECT   uvw[3];
440                                          /* record access */
# Line 374 | Line 443 | sumambient(            /* get interpolated ambient value */
443                  /*
444                   *  Ambient level test
445                   */
446 <                if (av->lvl > al)       /* list sorted, so this works */
446 >                if (av->lvl > al ||     /* list sorted, so this works */
447 >                                (av->lvl == al) & (av->weight < 0.9*r->rweight))
448                          break;
379 –                if (av->weight < 0.9*r->rweight)
380 –                        continue;
449                  /*
450                   *  Direction test using unperturbed normal
451                   */
# Line 389 | Line 457 | sumambient(            /* get interpolated ambient value */
457                  if (delta_r2 >= maxangle*maxangle)
458                          continue;
459                  /*
460 +                 *  Modified ray behind test
461 +                 */
462 +                VSUB(ck0, r->rop, av->pos);
463 +                d = DOT(ck0, uvw[2]);
464 +                if (d < -minarad*ambacc-.001)
465 +                        continue;
466 +                d /= av->rad[0];
467 +                delta_t2 = d*d;
468 +                if (delta_t2 >= ambacc*ambacc)
469 +                        continue;
470 +                /*
471                   *  Elliptical radii test based on Hessian
472                   */
473                  decodedir(uvw[0], av->udir);
474                  VCROSS(uvw[1], uvw[2], uvw[0]);
475 <                VSUB(ck0, av->pos, r->rop);
397 <                d = DOT(ck0, uvw[0]) / av->rad[0];
398 <                delta_t2 = d*d;
399 <                d = DOT(ck0, uvw[1]) / av->rad[1];
475 >                d = (u = DOT(ck0, uvw[0])) / av->rad[0];
476                  delta_t2 += d*d;
477 +                d = (v = DOT(ck0, uvw[1])) / av->rad[1];
478 +                delta_t2 += d*d;
479                  if (delta_t2 >= ambacc*ambacc)
480                          continue;
481                  /*
482 <                 *  Intersection behind test
482 >                 *  Test for potential light leak
483                   */
484 <                d = 0.0;
407 <                for (j = 0; j < 3; j++)
408 <                        d += (r->rop[j] - av->pos[j])*(uvw[2][j] + r->ron[j]);
409 <                if (d*0.5 < -minarad*ambacc-.001)
484 >                if (av->corral && plugaleak(r, av, uvw[2], atan2a(v,u)))
485                          continue;
486                  /*
487                   *  Extrapolate value and compute final weight (hat function)
# Line 418 | Line 493 | sumambient(            /* get interpolated ambient value */
493                  addcolor(acol, ct);
494                  wsum += d;
495          }
421 –        if (at->kid == NULL)
422 –                return(wsum);
423 –                                        /* sum children */
424 –        s *= 0.5;
425 –        for (i = 0; i < 8; i++) {
426 –                for (j = 0; j < 3; j++) {
427 –                        ck0[j] = c0[j];
428 –                        if (1<<j & i)
429 –                                ck0[j] += s;
430 –                        if (r->rop[j] < ck0[j] - OCTSCALE*s)
431 –                                break;
432 –                        if (r->rop[j] > ck0[j] + (1.0+OCTSCALE)*s)
433 –                                break;
434 –                }
435 –                if (j == 3)
436 –                        wsum += sumambient(acol, r, rn, al,
437 –                                                at->kid+i, ck0, s);
438 –        }
496          return(wsum);
497   }
498  
499  
500 < int
500 > static int
501   makeambient(            /* make a new ambient value for storage */
502          COLOR  acol,
503          RAY  *r,
# Line 459 | Line 516 | makeambient(           /* make a new ambient value for storage
516                  amb.weight = 1.25*r->rweight;
517          setcolor(acol, AVGREFL, AVGREFL, AVGREFL);
518                                                  /* compute ambient */
519 <        i = doambient(acol, r, amb.weight, uvw, amb.rad, amb.gpos, amb.gdir);
519 >        i = doambient(acol, r, amb.weight,
520 >                        uvw, amb.rad, amb.gpos, amb.gdir, &amb.corral);
521          scalecolor(acol, 1./AVGREFL);           /* undo assumed reflectance */
522 <        if (i <= 0)                             /* no Hessian => no storage */
522 >        if (i <= 0 || amb.rad[0] <= FTINY)      /* no Hessian or zero radius */
523                  return(i);
524                                                  /* store value */
525          VCOPY(amb.pos, r->rop);
# Line 479 | Line 537 | makeambient(           /* make a new ambient value for storage
537   }
538  
539  
540 < void
540 > static void
541   extambient(             /* extrapolate value at pv, nv */
542          COLOR  cr,
543          AMBVAL   *ap,
# Line 534 | Line 592 | avinsert(                              /* insert ambient value in our tree */
592          at = &atrunk;
593          VCOPY(ck0, thescene.cuorg);
594          s = thescene.cusize;
595 <        while (s*(OCTSCALE/2) > av->rad[1]*qambacc) {
595 >        while (s*(OCTSCALE/2) > av->rad[1]*ambacc) {
596                  if (at->kid == NULL)
597                          if ((at->kid = newambtree()) == NULL)
598                                  error(SYSTEM, "out of memory in avinsert");
# Line 549 | Line 607 | avinsert(                              /* insert ambient value in our tree */
607          }
608          avh.next = at->alist;           /* order by increasing level */
609          for (ap = &avh; ap->next != NULL; ap = ap->next)
610 <                if (ap->next->lvl >= av->lvl)
610 >                if ( ap->next->lvl > av->lvl ||
611 >                                (ap->next->lvl == av->lvl) &
612 >                                (ap->next->weight <= av->weight) )
613                          break;
614          av->next = ap->next;
615          ap->next = (AMBVAL*)av;
# Line 662 | Line 722 | sumambient(    /* get interpolated ambient value */
722                  /*
723                   *  Ambient level test.
724                   */
725 <                if (av->lvl > al)       /* list sorted, so this works */
725 >                if (av->lvl > al ||     /* list sorted, so this works */
726 >                                (av->lvl == al) & (av->weight < 0.9*r->rweight))
727                          break;
667 –                if (av->weight < 0.9*r->rweight)
668 –                        continue;
728                  /*
729                   *  Ambient radius test.
730                   */
# Line 847 | Line 906 | avinsert(                              /* insert ambient value in our tree */
906          }
907          avh.next = at->alist;           /* order by increasing level */
908          for (ap = &avh; ap->next != NULL; ap = ap->next)
909 <                if (ap->next->lvl >= av->lvl)
909 >                if ( ap->next->lvl > av->lvl ||
910 >                                (ap->next->lvl == av->lvl) &
911 >                                (ap->next->weight <= av->weight) )
912                          break;
913          av->next = ap->next;
914          ap->next = (AMBVAL*)av;
# Line 899 | Line 960 | avsave(                                /* insert and save an ambient value */
960          AMBVAL  *av
961   )
962   {
963 <        avinsert(avstore(av));
963 >        avstore(av);
964          if (ambfp == NULL)
965                  return;
966          if (writambval(av, ambfp) < 0)
# Line 914 | Line 975 | writerr:
975  
976  
977   static AMBVAL *
978 < avstore(                                /* allocate memory and store aval */
978 > avstore(                                /* allocate memory and save aval */
979          AMBVAL  *aval
980   )
981   {
# Line 932 | Line 993 | avstore(                               /* allocate memory and store aval */
993                  avsum += log(d);
994                  navsum++;
995          }
996 +        avinsert(av);                   /* insert in our cache tree */
997          return(av);
998   }
999  
# Line 1132 | Line 1194 | sortambvals(                   /* resort ambient values */
1194                  if (i_avlist < nambvals)
1195                          error(CONSISTENCY, "missing ambient values in sortambvals");
1196   #endif
1197 <                qsort((char *)avlist1, nambvals, sizeof(struct avl), &alatcmp);
1198 <                qsort((char *)avlist2, nambvals, sizeof(AMBVAL *), &aposcmp);
1197 >                qsort((char *)avlist1, nambvals, sizeof(struct avl), alatcmp);
1198 >                qsort((char *)avlist2, nambvals, sizeof(AMBVAL *), aposcmp);
1199                  for (i = 0; i < nambvals; i++) {
1200                          if (avlist1[i].p == NULL)
1201                                  continue;
# Line 1214 | Line 1276 | ambsync(void)                  /* synchronize ambient file */
1276                                  error(WARNING, errmsg);
1277                                  break;
1278                          }
1279 <                        avinsert(avstore(&avs));
1279 >                        avstore(&avs);
1280                          n -= AMBVALSIZ;
1281                  }
1282                  lastpos = flen - n;

Diff Legend

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