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

Comparing ray/src/rt/pmap.c (file contents):
Revision 2.14 by rschregle, Tue Mar 20 19:55:33 2018 UTC vs.
Revision 2.18 by rschregle, Fri Feb 19 02:10:35 2021 UTC

# Line 9 | Line 9 | static const char RCSid[] = "$Id$";
9  
10     Roland Schregle (roland.schregle@{hslu.ch, gmail.com})
11     (c) Fraunhofer Institute for Solar Energy Systems,
12 +       supported by the German Research Foundation
13 +       (DFG LU-204/10-2, "Fassadenintegrierte Regelsysteme FARESYS")
14     (c) Lucerne University of Applied Sciences and Arts,
15 <       supported by the Swiss National Science Foundation (SNSF, #147053)
15 >       supported by the Swiss National Science Foundation
16 >       (SNSF #147053, "Daylight Redirecting Components")
17     ======================================================================
18    
19     $Id$
# Line 25 | Line 28 | static const char RCSid[] = "$Id$";
28   #include "pmapbias.h"
29   #include "pmapdiag.h"
30   #include "otypes.h"
31 + #include "otspecial.h"
32   #include <time.h>
33   #if NIX
34     #include <sys/stat.h>
# Line 50 | Line 54 | static int photonParticipate (RAY *ray)
54     or 0 if absorbed and $*%&ed. Analogon to rayparticipate(). */
55   {
56     int i;
57 <   RREAL cosTheta, cosPhi, du, dv;
57 >   RREAL xi1, cosTheta, phi, du, dv;
58     const float cext = colorAvg(ray -> cext),
59 <               albedo = colorAvg(ray -> albedo);
59 >               albedo = colorAvg(ray -> albedo),
60 >               gecc = ray -> gecc, gecc2 = sqr(gecc);
61     FVECT u, v;
62     COLOR cvext;
63  
# Line 60 | Line 65 | static int photonParticipate (RAY *ray)
65     ray -> rmax = -log(pmapRandom(mediumState)) / cext;
66    
67     while (!localhit(ray, &thescene)) {
68 +      if (!incube(&thescene, ray -> rop)) {
69 +         /* Terminate photon if it has leaked from the scene */
70 + #ifdef DEBUG_PMAP
71 +         fprintf(stderr,
72 +                 "Volume photon leaked from scene at [%.3f %.3f %.3f]\n",
73 +                 ray -> rop [0], ray -> rop [1], ray -> rop [2]);
74 + #endif                
75 +         return 0;
76 +      }
77 +        
78        setcolor(cvext, exp(-ray -> rmax * ray -> cext [0]),
79                        exp(-ray -> rmax * ray -> cext [1]),
80                        exp(-ray -> rmax * ray -> cext [2]));
# Line 69 | Line 84 | static int photonParticipate (RAY *ray)
84        colorNorm(ray -> rcol);
85        VCOPY(ray -> rorg, ray -> rop);
86        
87 + #if 0
88        if (albedo > FTINY && ray -> rlvl > 0)
89 + #else
90 +      /* Store volume photons unconditionally in mist to also account for
91 +         direct inscattering from sources */
92 +      if (albedo > FTINY)
93 + #endif
94           /* Add to volume photon map */
95           newPhoton(volumePmap, ray);
96          
# Line 82 | Line 103 | static int photonParticipate (RAY *ray)
103        scalecolor(ray -> rcol, 1 / albedo);    
104        
105        /* Scatter photon */
106 <      cosTheta = ray -> gecc <= FTINY ? 2 * pmapRandom(scatterState) - 1
107 <                                      : 1 / (2 * ray -> gecc) *
108 <                                            (1 + ray -> gecc * ray -> gecc -
109 <                                               (1 - ray -> gecc * ray -> gecc) /
110 <                                               (1 - ray -> gecc + 2 * ray -> gecc *
111 <                                                  pmapRandom(scatterState)));
112 <                                                  
113 <      cosPhi = cos(2 * PI * pmapRandom(scatterState));
114 <      du = dv = sqrt(1 - cosTheta * cosTheta);   /* sin(theta) */
115 <      du *= cosPhi;
116 <      dv *= sqrt(1 - cosPhi * cosPhi);           /* sin(phi) */
106 >      xi1 = pmapRandom(scatterState);
107 >      cosTheta = ray -> gecc <= FTINY
108 >                    ? 2 * xi1 - 1
109 >                    : 0.5 / gecc *
110 >                      (1 + gecc2 - sqr((1 - gecc2) /
111 >                                       (1 + gecc * (2 * xi1 - 1))));
112 >
113 >      phi = 2 * PI * pmapRandom(scatterState);
114 >      du = dv = sqrt(1 - sqr(cosTheta));     /* sin(theta) */
115 >      du *= cos(phi);
116 >      dv *= sin(phi);
117        
118        /* Get axes u & v perpendicular to photon direction */
119        i = 0;
# Line 106 | Line 127 | static int photonParticipate (RAY *ray)
127        for (i = 0; i < 3; i++)
128           ray -> rdir [i] = du * u [i] + dv * v [i] +
129                             cosTheta * ray -> rdir [i];
130 +
131        ray -> rlvl++;
132        ray -> rmax = -log(pmapRandom(mediumState)) / cext;
133     }  
134 <    
134 >  
135 >   /* Passed through medium until intersecting local object */  
136     setcolor(cvext, exp(-ray -> rot * ray -> cext [0]),
137                     exp(-ray -> rot * ray -> cext [1]),
138                     exp(-ray -> rot * ray -> cext [2]));
139                    
140     /* Modify ray color and normalise */
141     multcolor(ray -> rcol, cvext);
142 <   colorNorm(ray -> rcol);
143 <  
121 <   /* Passed through medium */  
142 >   colorNorm(ray -> rcol);  
143 >
144     return 1;
145   }
146  
# Line 151 | Line 173 | void tracePhoton (RAY *ray)
173     if (localhit(ray, &thescene)) {
174        mod = ray -> ro -> omod;
175  
176 <      if (port && ray -> ro != port) {
176 >      /* XXX: Is port -> omod != mod sufficient here? Probably not... */
177 >      if (
178 >         port && ray -> ro != port &&
179 >         findmaterial(port) != findmaterial(ray -> ro)
180 >      ) {
181           /* !!! PHOTON PORT REJECTION SAMPLING HACK !!!
182 <          * Terminate photon if emitted from port without intersecting it;
183 <          * this can happen when the port's partitions extend beyond its
182 >          * Terminate photon if emitted from port without intersecting it or
183 >          * its other associated surfaces or same material.
184 >          * This can happen when the port's partitions extend beyond its
185            * actual geometry, e.g.  with polygons.  Since the total flux
186            * relayed by the port is based on the (in this case) larger
187            * partition area, it is overestimated; terminating these photons
# Line 444 | Line 471 | void distribPhotons (PhotonMap **pmaps, unsigned numPr
471           pmapSeed(randSeed + (proc + 3) % numProc, mediumState);
472           pmapSeed(randSeed + (proc + 4) % numProc, scatterState);
473           pmapSeed(randSeed + (proc + 5) % numProc, rouletteState);
474 <                  
474 >              
475 > #ifdef DEBUG_PMAP          
476 >         /* Output child process PID after random delay to prevent corrupted
477 >          * console output due to race condition */
478 >         usleep(1e6 * pmapRandom(rouletteState));
479 >         fprintf(stderr, "Proc %d: PID = %d "
480 >                 "(waiting 10 sec to attach debugger...)\n",
481 >                 proc, getpid());
482 >         /* Allow time for debugger to attach to child process */
483 >         sleep(10);
484 > #endif            
485 >
486           for (t = 0; t < NUM_PMAP_TYPES; t++)
487              lastNumPhotons [t] = 0;
488              

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines