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

Comparing ray/src/rt/pmapdata.c (file contents):
Revision 2.1 by greg, Tue Feb 24 19:39:26 2015 UTC vs.
Revision 2.14 by rschregle, Thu Feb 4 11:36:59 2016 UTC

# Line 1 | Line 1
1 + #ifndef lint
2 + static const char RCSid[] = "$Id$";
3 + #endif
4   /*
5     ==================================================================
6     Photon map data structures and kd-tree handling
7  
8     Roland Schregle (roland.schregle@{hslu.ch, gmail.com})
9     (c) Fraunhofer Institute for Solar Energy Systems,
10 <       Lucerne University of Applied Sciences & Arts
10 >   (c) Lucerne University of Applied Sciences and Arts,
11 >   supported by the Swiss National Science Foundation (SNSF, #147053)
12     ==================================================================  
13    
10   $Id$
14   */
15  
16  
# Line 18 | Line 21
21   #include "otypes.h"
22   #include "source.h"
23   #include "rcontrib.h"
24 + #include "random.h"
25  
26  
27  
# Line 65 | Line 69 | void initPhotonMap (PhotonMap *pmap, PhotonMapType t)
69   const PhotonPrimary* addPhotonPrimary (PhotonMap *pmap, const RAY *ray)
70   {
71     PhotonPrimary *prim = NULL;
72 +   FVECT dvec;
73    
74     if (!pmap || !ray)
75        return NULL;
# Line 97 | Line 102 | const PhotonPrimary* addPhotonPrimary (PhotonMap *pmap
102     prim -> srcIdx = -1;
103        
104     /* Reverse incident direction to point to light source */
105 <   prim -> dir [0] = -ray -> rdir [0];
106 <   prim -> dir [1] = -ray -> rdir [1];
107 <   prim -> dir [2] = -ray -> rdir [2];
105 >   dvec [0] = -ray -> rdir [0];
106 >   dvec [1] = -ray -> rdir [1];
107 >   dvec [2] = -ray -> rdir [2];
108 >   prim -> dir = encodedir(dvec);
109  
110 <   VCOPY(prim -> org, ray -> rorg);
110 >   VCOPY(prim -> pos, ray -> rop);
111    
112     return prim;
113   }
# Line 225 | Line 231 | static void nearestNeighbours (PhotonMap* pmap, const
231           nearestNeighbours(pmap, pos, norm, node << 1);
232     }
233  
234 <   /* Reject photon if normal faces away (ignored for volume photons) */
235 <   if (norm && DOT(norm, p -> norm) <= 0)
234 >   /* Reject photon if normal faces away (ignored for volume photons) with
235 >    * 50% tolerance to account for perturbation; note photon normal is coded
236 >    * in range [-127,127].  */
237 >   if (norm && DOT(norm, p -> norm) <= 63.5 * frandom())
238        return;
239        
240     if (isContribPmap(pmap) && pmap -> srcContrib) {
# Line 237 | Line 245 | static void nearestNeighbours (PhotonMap* pmap, const
245        if (srcIdx < 0 || srcIdx >= nsources)
246           error(INTERNAL, "invalid light source index in photon map");
247        
248 <      srcMod = objptr(source [srcIdx].so -> omod);
248 >      srcMod = findmaterial(source [srcIdx].so);
249  
250        /* Reject photon if contributions from light source which emitted it
251         * are not sought */
# Line 308 | Line 316 | static void nearestNeighbours (PhotonMap* pmap, const
316   /* Threshold below which we assume increasing max radius won't help */
317   #define PMAP_SHORT_LOOKUP_THRESH 1
318  
319 + /* Coefficient for adaptive maximum search radius */
320 + #define PMAP_MAXDIST_COEFF 100
321 +
322 +
323   void findPhotons (PhotonMap* pmap, const RAY* ray)
324   {
325     float pos [3], norm [3];
# Line 329 | Line 341 | void findPhotons (PhotonMap* pmap, const RAY* ray)
341        pmap -> minError = FHUGE;
342        pmap -> maxError = -FHUGE;
343        pmap -> rmsError = 0;
344 <      /* Maximum search radius limit based on avg photon distance to
345 <       * centre of gravity */
344 >      /* SQUARED max search radius limit is based on avg photon distance to
345 >       * centre of gravity, unless fixed by user (maxDistFix > 0) */
346        pmap -> maxDist0 = pmap -> maxDistLimit =
347 <         maxDistCoeff * pmap -> squeueSize * pmap -> CoGdist /
348 <         pmap -> heapSize;
347 >         maxDistFix > 0 ? maxDistFix * maxDistFix
348 >                        : PMAP_MAXDIST_COEFF * pmap -> squeueSize *
349 >                          pmap -> CoGdist / pmap -> heapSize;
350     }
351  
352     do {
# Line 352 | Line 365 | void findPhotons (PhotonMap* pmap, const RAY* ray)
365           VCOPY(norm, ray -> ron);
366           nearestNeighbours(pmap, pos, norm, 1);
367        }
368 <      
368 >
369 > #ifdef PMAP_ITSYBITSY
370 >      if (pmap -> maxDist < FTINY) {
371 >         sprintf(errmsg, "itsy bitsy teeny weeny photon search radius %e",
372 >                 sqrt(pmap -> maxDist));
373 >         error(WARNING, errmsg);
374 >      }
375 > #endif
376 >
377        if (pmap -> squeueEnd < pmap -> squeueSize * pmap -> gatherTolerance) {
378           /* Short lookup; too few photons found */
379           if (pmap -> squeueEnd > PMAP_SHORT_LOOKUP_THRESH) {
# Line 368 | Line 389 | void findPhotons (PhotonMap* pmap, const RAY* ray)
389                      ray -> ro ? ray -> ro -> oname : "<null>");
390              error(WARNING, errmsg);        
391   #endif            
392 <            
392 >
393 >            /* Bail out after warning if maxDist is fixed */
394 >            if (maxDistFix > 0)
395 >               return;
396 >            
397              if (pmap -> maxDist0 < pmap -> maxDistLimit) {
398                 /* Increase max search radius if below limit & redo search */
399                 pmap -> maxDist0 *= PMAP_MAXDIST_INC;
400   #ifdef PMAP_LOOKUP_REDO
401                 redo = 1;
402   #endif              
378
403   #ifdef PMAP_LOOKUP_WARN
404                 sprintf(errmsg,
405                         redo ? "restarting photon lookup with max radius %.1e"
406                              : "max photon lookup radius adjusted to %.1e",
407 <                       pmap -> maxDist0);
407 >                       sqrt(pmap -> maxDist0));
408                 error(WARNING, errmsg);
409   #endif
410              }
411 + #ifdef PMAP_LOOKUP_REDO
412              else {
413                 sprintf(errmsg, "max photon lookup radius clamped to %.1e",
414 <                       pmap -> maxDist0);
414 >                       sqrt(pmap -> maxDist0));
415                 error(WARNING, errmsg);
416              }
417 + #endif
418           }
419          
420           /* Reset successful lookup counter */
421           pmap -> numLookups = 0;
422 <      }  
422 >      }
423        else {
424 +         /* Bail out after warning if maxDist is fixed */
425 +         if (maxDistFix > 0)
426 +            return;
427 +            
428           /* Increment successful lookup counter and reduce max search radius if
429            * wraparound */
430           pmap -> numLookups = (pmap -> numLookups + 1) % PMAP_MAXDIST_CNT;
# Line 403 | Line 433 | void findPhotons (PhotonMap* pmap, const RAY* ray)
433              
434           redo = 0;
435        }
436 +      
437     } while (redo);
438   }
439  
# Line 442 | Line 473 | static void nearest1Neighbour (PhotonMap *pmap, const
473     dv [2] = pos [2] - p -> pos [2];
474     d2 = DOT(dv, dv);
475    
476 <   if (d2 < pmap -> maxDist && DOT(norm, p -> norm) > 0) {
477 <      /* Closest photon so far with similar normal */
476 >   if (d2 < pmap -> maxDist && DOT(norm, p -> norm) > 63.5 * frandom()) {
477 >      /* Closest photon so far with similar normal. We allow a 50% tolerance
478 >       * to account for perturbation in the latter; note the photon normal
479 >       * is coded in the range [-127,127].  */  
480        pmap -> maxDist = d2;
481        *photon = p;
482     }

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines