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.11 by rschregle, Tue May 17 17:39:47 2016 UTC vs.
Revision 2.12 by greg, Mon Sep 26 20:19:30 2016 UTC

# Line 30 | Line 30 | static const char RCSid[] = "$Id$";
30   #include <sys/mman.h>
31   #include <sys/wait.h>
32  
33 – #define PMAP_REV  "$Revision$"
33  
34  
36 – extern char *octname;
37 –
38 –
39 –
40 – /* Photon map lookup functions per type */
41 – void (*pmapLookup [NUM_PMAP_TYPES])(PhotonMap*, RAY*, COLOR) = {
42 –   photonDensity, photonPreCompDensity, photonDensity, volumePhotonDensity,
43 –   photonDensity, NULL
44 – };
45 –
46 –
47 –
48 – void colorNorm (COLOR c)
49 – /* Normalise colour channels to average of 1 */
50 – {
51 –   const float avg = colorAvg(c);
52 –  
53 –   if (!avg)
54 –      return;
55 –      
56 –   c [0] /= avg;
57 –   c [1] /= avg;
58 –   c [2] /= avg;
59 – }
60 –
61 –
62 –
63 – void loadPmaps (PhotonMap **pmaps, const PhotonMapParams *parm)
64 – {
65 –   unsigned t;
66 –   struct stat octstat, pmstat;
67 –   PhotonMap *pm;
68 –   PhotonMapType type;
69 –  
70 –   for (t = 0; t < NUM_PMAP_TYPES; t++)
71 –      if (setPmapParam(&pm, parm + t)) {        
72 –         /* Check if photon map newer than octree */
73 –         if (pm -> fileName && octname &&
74 –             !stat(pm -> fileName, &pmstat) && !stat(octname, &octstat) &&
75 –             octstat.st_mtime > pmstat.st_mtime) {
76 –            sprintf(errmsg, "photon map in file %s may be stale",
77 –                    pm -> fileName);
78 –            error(USER, errmsg);
79 –         }
80 –        
81 –         /* Load photon map from file and get its type */
82 –         if ((type = loadPhotonMap(pm, pm -> fileName)) == PMAP_TYPE_NONE)
83 –            error(USER, "failed loading photon map");
84 –            
85 –         /* Assign to appropriate photon map type (deleting previously
86 –          * loaded photon map of same type if necessary) */
87 –         if (pmaps [type]) {
88 –            deletePhotons(pmaps [type]);
89 –            free(pmaps [type]);
90 –         }
91 –         pmaps [type] = pm;
92 –        
93 –         /* Check for invalid density estimate bandwidth */                            
94 –         if (pm -> maxGather > pm -> numPhotons) {
95 –            error(WARNING, "adjusting density estimate bandwidth");
96 –            pm -> minGather = pm -> maxGather = pm -> numPhotons;
97 –         }
98 –      }
99 – }
100 –
101 –
102 –
35   void savePmaps (const PhotonMap **pmaps, int argc, char **argv)
36   {
37     unsigned t;
# Line 111 | Line 43 | void savePmaps (const PhotonMap **pmaps, int argc, cha
43   }                  
44  
45  
114 –  
115 – void cleanUpPmaps (PhotonMap **pmaps)
116 – {
117 –   unsigned t;
118 –  
119 –   for (t = 0; t < NUM_PMAP_TYPES; t++) {
120 –      if (pmaps [t]) {
121 –         deletePhotons(pmaps [t]);
122 –         free(pmaps [t]);
123 –      }
124 –   }
125 – }
126 –
127 –
46      
47   static int photonParticipate (RAY *ray)
48   /* Trace photon through participating medium. Returns 1 if passed through,
# Line 778 | Line 696 | void distribPhotons (PhotonMap **pmaps, unsigned numPr
696     /* Precompute photon irradiance if necessary */
697     if (preCompPmap)
698        preComputeGlobal(preCompPmap);
781 – }
782 –
783 –
784 –
785 – void photonDensity (PhotonMap *pmap, RAY *ray, COLOR irrad)
786 – /* Photon density estimate. Returns irradiance at ray -> rop. */
787 – {
788 –   unsigned                      i;
789 –   float                         r;
790 –   COLOR                         flux;
791 –   Photon                        *photon;
792 –   const PhotonSearchQueueNode   *sqn;
793 –
794 –   setcolor(irrad, 0, 0, 0);
795 –
796 –   if (!pmap -> maxGather)
797 –      return;
798 –      
799 –   /* Ignore sources */
800 –   if (ray -> ro && islight(objptr(ray -> ro -> omod) -> otype))
801 –      return;
802 –        
803 –   findPhotons(pmap, ray);
804 –  
805 –   /* Need at least 2 photons */
806 –   if (pmap -> squeue.tail < 2) {
807 – #ifdef PMAP_NONEFOUND  
808 –      sprintf(errmsg, "no photons found on %s at (%.3f, %.3f, %.3f)",
809 –              ray -> ro ? ray -> ro -> oname : "<null>",
810 –              ray -> rop [0], ray -> rop [1], ray -> rop [2]);
811 –      error(WARNING, errmsg);
812 – #endif      
813 –
814 –      return;
815 –   }
816 –
817 –   if (pmap -> minGather == pmap -> maxGather) {
818 –      /* No bias compensation. Just do a plain vanilla estimate */
819 –      sqn = pmap -> squeue.node + 1;
820 –      
821 –      /* Average radius between furthest two photons to improve accuracy */      
822 –      r = max(sqn -> dist2, (sqn + 1) -> dist2);
823 –      r = 0.25 * (pmap -> maxDist2 + r + 2 * sqrt(pmap -> maxDist2 * r));  
824 –      
825 –      /* Skip the extra photon */
826 –      for (i = 1 ; i < pmap -> squeue.tail; i++, sqn++) {
827 –         photon = getNearestPhoton(&pmap -> squeue, sqn -> idx);
828 –         getPhotonFlux(photon, flux);        
829 – #ifdef PMAP_EPANECHNIKOV
830 –         /* Apply Epanechnikov kernel to photon flux based on photon dist */
831 –         scalecolor(flux, 2 * (1 - sqn -> dist2 / r));
832 – #endif  
833 –         addcolor(irrad, flux);
834 –      }
835 –      
836 –      /* Divide by search area PI * r^2, 1 / PI required as ambient
837 –         normalisation factor */        
838 –      scalecolor(irrad, 1 / (PI * PI * r));
839 –      
840 –      return;
841 –   }
842 –   else
843 –      /* Apply bias compensation to density estimate */
844 –      biasComp(pmap, irrad);
845 – }
846 –
847 –
848 –
849 – void photonPreCompDensity (PhotonMap *pmap, RAY *r, COLOR irrad)
850 – /* Returns precomputed photon density estimate at ray -> rop. */
851 – {
852 –   Photon p;
853 –  
854 –   setcolor(irrad, 0, 0, 0);
855 –
856 –   /* Ignore sources */
857 –   if (r -> ro && islight(objptr(r -> ro -> omod) -> otype))
858 –      return;
859 –      
860 –   find1Photon(preCompPmap, r, &p);
861 –   getPhotonFlux(&p, irrad);
862 – }
863 –
864 –
865 –
866 – void volumePhotonDensity (PhotonMap *pmap, RAY *ray, COLOR irrad)
867 – /* Photon volume density estimate. Returns irradiance at ray -> rop. */
868 – {
869 –   unsigned                      i;
870 –   float                         r, gecc2, ph;
871 –   COLOR                         flux;
872 –   Photon                        *photon;
873 –   const PhotonSearchQueueNode   *sqn;
874 –
875 –   setcolor(irrad, 0, 0, 0);
876 –  
877 –   if (!pmap -> maxGather)
878 –      return;
879 –      
880 –   findPhotons(pmap, ray);
881 –  
882 –   /* Need at least 2 photons */
883 –   if (pmap -> squeue.tail < 2)
884 –      return;
885 –
886 – #if 0      
887 –   /* Volume biascomp disabled (probably redundant) */
888 –   if (pmap -> minGather == pmap -> maxGather)
889 – #endif  
890 –   {
891 –      /* No bias compensation. Just do a plain vanilla estimate */
892 –      gecc2 = ray -> gecc * ray -> gecc;
893 –      sqn = pmap -> squeue.node + 1;
894 –      
895 –      /* Average radius between furthest two photons to improve accuracy */      
896 –      r = max(sqn -> dist2, (sqn + 1) -> dist2);
897 –      r = 0.25 * (pmap -> maxDist2 + r + 2 * sqrt(pmap -> maxDist2 * r));  
898 –      
899 –      /* Skip the extra photon */
900 –      for (i = 1; i < pmap -> squeue.tail; i++, sqn++) {
901 –         photon = getNearestPhoton(&pmap -> squeue, sqn -> idx);
902 –        
903 –         /* Compute phase function for inscattering from photon */
904 –         if (gecc2 <= FTINY)
905 –            ph = 1;
906 –         else {
907 –            ph = DOT(ray -> rdir, photon -> norm) / 127;
908 –            ph = 1 + gecc2 - 2 * ray -> gecc * ph;
909 –            ph = (1 - gecc2) / (ph * sqrt(ph));
910 –         }
911 –        
912 –         getPhotonFlux(photon, flux);
913 –         scalecolor(flux, ph);
914 –         addcolor(irrad, flux);
915 –      }
916 –      
917 –      /* Divide by search volume 4 / 3 * PI * r^3 and phase function
918 –         normalization factor 1 / (4 * PI) */
919 –      scalecolor(irrad, 3 / (16 * PI * PI * r * sqrt(r)));
920 –      return;
921 –   }
922 – #if 0
923 –   else
924 –      /* Apply bias compensation to density estimate */
925 –      volumeBiasComp(pmap, ray, irrad);
926 – #endif      
699   }

Diff Legend

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