--- ray/src/rt/pmapdata.c	2015/04/22 15:44:57	2.2
+++ ray/src/rt/pmapdata.c	2015/05/26 13:31:19	2.8
@@ -4,10 +4,11 @@
 
    Roland Schregle (roland.schregle@{hslu.ch, gmail.com})
    (c) Fraunhofer Institute for Solar Energy Systems,
-       Lucerne University of Applied Sciences & Arts
+   (c) Lucerne University of Applied Sciences and Arts,
+   supported by the Swiss National Science Foundation (SNSF, #147053)
    ==================================================================   
    
-   $Id: pmapdata.c,v 2.2 2015/04/22 15:44:57 rschregle Exp $
+   $Id: pmapdata.c,v 2.8 2015/05/26 13:31:19 rschregle Exp $
 */
 
 
@@ -18,6 +19,7 @@
 #include "otypes.h"
 #include "source.h"
 #include "rcontrib.h"
+#include "random.h"
 
 
 
@@ -65,6 +67,7 @@ void initPhotonMap (PhotonMap *pmap, PhotonMapType t)
 const PhotonPrimary* addPhotonPrimary (PhotonMap *pmap, const RAY *ray)
 {
    PhotonPrimary *prim = NULL;
+   FVECT dvec;
    
    if (!pmap || !ray)
       return NULL;
@@ -97,11 +100,12 @@ const PhotonPrimary* addPhotonPrimary (PhotonMap *pmap
    prim -> srcIdx = -1;
       
    /* Reverse incident direction to point to light source */
-   prim -> dir [0] = -ray -> rdir [0];
-   prim -> dir [1] = -ray -> rdir [1];
-   prim -> dir [2] = -ray -> rdir [2];
+   dvec [0] = -ray -> rdir [0];
+   dvec [1] = -ray -> rdir [1];
+   dvec [2] = -ray -> rdir [2];
+   prim -> dir = encodedir(dvec);
 
-   VCOPY(prim -> org, ray -> rorg);
+   VCOPY(prim -> pos, ray -> rop);
    
    return prim;
 }
@@ -226,7 +230,7 @@ static void nearestNeighbours (PhotonMap* pmap, const 
    }
 
    /* Reject photon if normal faces away (ignored for volume photons) */
-   if (norm && DOT(norm, p -> norm) <= 0)
+   if (norm && DOT(norm, p -> norm) <= 0.5 * frandom())
       return;
       
    if (isContribPmap(pmap) && pmap -> srcContrib) {
@@ -237,7 +241,7 @@ static void nearestNeighbours (PhotonMap* pmap, const 
       if (srcIdx < 0 || srcIdx >= nsources)
          error(INTERNAL, "invalid light source index in photon map");
       
-      srcMod = objptr(source [srcIdx].so -> omod);
+      srcMod = findmaterial(source [srcIdx].so);
 
       /* Reject photon if contributions from light source which emitted it
        * are not sought */
@@ -308,6 +312,10 @@ static void nearestNeighbours (PhotonMap* pmap, const 
 /* Threshold below which we assume increasing max radius won't help */
 #define PMAP_SHORT_LOOKUP_THRESH 1
 
+/* Coefficient for adaptive maximum search radius */
+#define PMAP_MAXDIST_COEFF 100
+
+
 void findPhotons (PhotonMap* pmap, const RAY* ray)
 {
    float pos [3], norm [3];
@@ -329,11 +337,12 @@ void findPhotons (PhotonMap* pmap, const RAY* ray)
       pmap -> minError = FHUGE;
       pmap -> maxError = -FHUGE;
       pmap -> rmsError = 0;
-      /* Maximum search radius limit based on avg photon distance to
-       * centre of gravity */
+      /* Maximum search radius limit is based on avg photon distance to
+       * centre of gravity, unless fixed by user (maxDistFix > 0) */
       pmap -> maxDist0 = pmap -> maxDistLimit = 
-         maxDistCoeff * pmap -> squeueSize * pmap -> CoGdist / 
-         pmap -> heapSize;
+         maxDistFix > 0 ? maxDistFix
+                        : PMAP_MAXDIST_COEFF * pmap -> squeueSize * 
+                          pmap -> CoGdist / pmap -> heapSize;
    }
 
    do {
@@ -352,7 +361,7 @@ void findPhotons (PhotonMap* pmap, const RAY* ray)
          VCOPY(norm, ray -> ron);
          nearestNeighbours(pmap, pos, norm, 1);
       }
-      
+
       if (pmap -> squeueEnd < pmap -> squeueSize * pmap -> gatherTolerance) {
          /* Short lookup; too few photons found */
          if (pmap -> squeueEnd > PMAP_SHORT_LOOKUP_THRESH) {
@@ -368,14 +377,17 @@ void findPhotons (PhotonMap* pmap, const RAY* ray)
                     ray -> ro ? ray -> ro -> oname : "<null>");
             error(WARNING, errmsg);         
 #endif             
-             
+
+            /* Bail out after warning if maxDist is fixed */
+            if (maxDistFix > 0)
+               return;
+            
             if (pmap -> maxDist0 < pmap -> maxDistLimit) {
                /* Increase max search radius if below limit & redo search */
                pmap -> maxDist0 *= PMAP_MAXDIST_INC;
 #ifdef PMAP_LOOKUP_REDO
                redo = 1;
 #endif               
-
 #ifdef PMAP_LOOKUP_WARN
                sprintf(errmsg, 
                        redo ? "restarting photon lookup with max radius %.1e"
@@ -395,8 +407,12 @@ void findPhotons (PhotonMap* pmap, const RAY* ray)
          
          /* Reset successful lookup counter */
          pmap -> numLookups = 0;
-      }  
+      }
       else {
+         /* Bail out after warning if maxDist is fixed */
+         if (maxDistFix > 0)
+            return;
+            
          /* Increment successful lookup counter and reduce max search radius if
           * wraparound */
          pmap -> numLookups = (pmap -> numLookups + 1) % PMAP_MAXDIST_CNT;
@@ -405,6 +421,7 @@ void findPhotons (PhotonMap* pmap, const RAY* ray)
             
          redo = 0;
       }
+      
    } while (redo);
 }
 
@@ -444,7 +461,7 @@ static void nearest1Neighbour (PhotonMap *pmap, const 
    dv [2] = pos [2] - p -> pos [2];
    d2 = DOT(dv, dv);
    
-   if (d2 < pmap -> maxDist && DOT(norm, p -> norm) > 0) {
+   if (d2 < pmap -> maxDist && DOT(norm, p -> norm) > 0.5 * frandom()) {
       /* Closest photon so far with similar normal */
       pmap -> maxDist = d2;
       *photon = p;