ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/rt/pmap.c
Revision: 2.15
Committed: Thu Jun 7 19:26:04 2018 UTC (6 years, 10 months ago) by rschregle
Content type: text/plain
Branch: MAIN
CVS Tags: rad5R2
Changes since 2.14: +20 -3 lines
Log Message:
Added photon distrib PID output to attach debugger

File Contents

# User Rev Content
1 greg 2.9 #ifndef lint
2 rschregle 2.15 static const char RCSid[] = "$Id: pmap.c,v 2.14 2018/03/20 19:55:33 rschregle Exp $";
3 greg 2.9 #endif
4 rschregle 2.11
5 rschregle 2.13
6 greg 2.1 /*
7 rschregle 2.11 ======================================================================
8 greg 2.1 Photon map main module
9    
10     Roland Schregle (roland.schregle@{hslu.ch, gmail.com})
11     (c) Fraunhofer Institute for Solar Energy Systems,
12 rschregle 2.4 (c) Lucerne University of Applied Sciences and Arts,
13 rschregle 2.11 supported by the Swiss National Science Foundation (SNSF, #147053)
14     ======================================================================
15 greg 2.1
16 rschregle 2.15 $Id: pmap.c,v 2.14 2018/03/20 19:55:33 rschregle Exp $
17 greg 2.1 */
18    
19    
20     #include "pmap.h"
21     #include "pmapmat.h"
22     #include "pmapsrc.h"
23     #include "pmaprand.h"
24     #include "pmapio.h"
25     #include "pmapbias.h"
26     #include "pmapdiag.h"
27     #include "otypes.h"
28     #include <time.h>
29 rschregle 2.13 #if NIX
30     #include <sys/stat.h>
31     #include <sys/mman.h>
32     #include <sys/wait.h>
33     #endif
34 greg 2.1
35    
36     void savePmaps (const PhotonMap **pmaps, int argc, char **argv)
37     {
38     unsigned t;
39    
40     for (t = 0; t < NUM_PMAP_TYPES; t++) {
41     if (pmaps [t])
42 greg 2.7 savePhotonMap(pmaps [t], pmaps [t] -> fileName, argc, argv);
43 greg 2.1 }
44     }
45    
46    
47    
48     static int photonParticipate (RAY *ray)
49     /* Trace photon through participating medium. Returns 1 if passed through,
50     or 0 if absorbed and $*%&ed. Analogon to rayparticipate(). */
51     {
52     int i;
53     RREAL cosTheta, cosPhi, du, dv;
54     const float cext = colorAvg(ray -> cext),
55     albedo = colorAvg(ray -> albedo);
56     FVECT u, v;
57     COLOR cvext;
58    
59     /* Mean free distance until interaction with medium */
60     ray -> rmax = -log(pmapRandom(mediumState)) / cext;
61    
62     while (!localhit(ray, &thescene)) {
63     setcolor(cvext, exp(-ray -> rmax * ray -> cext [0]),
64     exp(-ray -> rmax * ray -> cext [1]),
65     exp(-ray -> rmax * ray -> cext [2]));
66    
67     /* Modify ray color and normalise */
68     multcolor(ray -> rcol, cvext);
69     colorNorm(ray -> rcol);
70     VCOPY(ray -> rorg, ray -> rop);
71    
72 rschregle 2.15 #if 0
73 rschregle 2.11 if (albedo > FTINY && ray -> rlvl > 0)
74 rschregle 2.15 #else
75     /* Store volume photons unconditionally in mist to also account for
76     direct inscattering from sources */
77     if (albedo > FTINY)
78     #endif
79 greg 2.1 /* Add to volume photon map */
80 rschregle 2.11 newPhoton(volumePmap, ray);
81 greg 2.1
82     /* Absorbed? */
83 rschregle 2.11 if (pmapRandom(rouletteState) > albedo)
84     return 0;
85 greg 2.1
86     /* Colour bleeding without attenuation (?) */
87     multcolor(ray -> rcol, ray -> albedo);
88     scalecolor(ray -> rcol, 1 / albedo);
89    
90     /* Scatter photon */
91     cosTheta = ray -> gecc <= FTINY ? 2 * pmapRandom(scatterState) - 1
92     : 1 / (2 * ray -> gecc) *
93     (1 + ray -> gecc * ray -> gecc -
94     (1 - ray -> gecc * ray -> gecc) /
95     (1 - ray -> gecc + 2 * ray -> gecc *
96     pmapRandom(scatterState)));
97    
98     cosPhi = cos(2 * PI * pmapRandom(scatterState));
99     du = dv = sqrt(1 - cosTheta * cosTheta); /* sin(theta) */
100     du *= cosPhi;
101     dv *= sqrt(1 - cosPhi * cosPhi); /* sin(phi) */
102    
103     /* Get axes u & v perpendicular to photon direction */
104     i = 0;
105     do {
106     v [0] = v [1] = v [2] = 0;
107     v [i++] = 1;
108     fcross(u, v, ray -> rdir);
109     } while (normalize(u) < FTINY);
110     fcross(v, ray -> rdir, u);
111    
112     for (i = 0; i < 3; i++)
113     ray -> rdir [i] = du * u [i] + dv * v [i] +
114     cosTheta * ray -> rdir [i];
115     ray -> rlvl++;
116     ray -> rmax = -log(pmapRandom(mediumState)) / cext;
117     }
118    
119     setcolor(cvext, exp(-ray -> rot * ray -> cext [0]),
120     exp(-ray -> rot * ray -> cext [1]),
121     exp(-ray -> rot * ray -> cext [2]));
122    
123     /* Modify ray color and normalise */
124     multcolor(ray -> rcol, cvext);
125     colorNorm(ray -> rcol);
126    
127     /* Passed through medium */
128     return 1;
129     }
130    
131    
132    
133     void tracePhoton (RAY *ray)
134     /* Follow photon as it bounces around... */
135     {
136     long mod;
137 rschregle 2.13 OBJREC *mat, *port = NULL;
138    
139     if (!ray -> parent) {
140     /* !!! PHOTON PORT REJECTION SAMPLING HACK: get photon port for
141     * !!! primary ray from ray -> ro, then reset the latter to NULL so
142     * !!! as not to interfere with localhit() */
143     port = ray -> ro;
144     ray -> ro = NULL;
145     }
146 greg 2.1
147     if (ray -> rlvl > photonMaxBounce) {
148 rschregle 2.5 #ifdef PMAP_RUNAWAY_WARN
149 greg 2.1 error(WARNING, "runaway photon!");
150 rschregle 2.5 #endif
151 greg 2.1 return;
152     }
153 rschregle 2.5
154 greg 2.1 if (colorAvg(ray -> cext) > FTINY && !photonParticipate(ray))
155     return;
156 rschregle 2.13
157 greg 2.1 if (localhit(ray, &thescene)) {
158     mod = ray -> ro -> omod;
159 rschregle 2.13
160     if (port && ray -> ro != port) {
161     /* !!! PHOTON PORT REJECTION SAMPLING HACK !!!
162     * Terminate photon if emitted from port without intersecting it;
163     * this can happen when the port's partitions extend beyond its
164     * actual geometry, e.g. with polygons. Since the total flux
165     * relayed by the port is based on the (in this case) larger
166     * partition area, it is overestimated; terminating these photons
167     * constitutes rejection sampling and thereby compensates any bias
168     * incurred by the overestimated flux. */
169     #ifdef PMAP_PORTREJECT_WARN
170     sprintf(errmsg, "photon outside port %s", ray -> ro -> oname);
171     error(WARNING, errmsg);
172     #endif
173     return;
174     }
175    
176 greg 2.1 if ((ray -> clipset && inset(ray -> clipset, mod)) || mod == OVOID) {
177     /* Transfer ray if modifier is VOID or clipped within antimatta */
178     RAY tray;
179     photonRay(ray, &tray, PMAP_XFER, NULL);
180     tracePhoton(&tray);
181     }
182     else {
183     /* Scatter for modifier material */
184     mat = objptr(mod);
185     photonScatter [mat -> otype] (mat, ray);
186     }
187     }
188     }
189    
190    
191    
192     static void preComputeGlobal (PhotonMap *pmap)
193 rschregle 2.11 /* Precompute irradiance from global photons for final gathering for
194     a random subset of finalGather * pmap -> numPhotons photons, and builds
195     the photon map, discarding the original photons. */
196     /* !!! NOTE: PRECOMPUTATION WITH OOC CURRENTLY WITHOUT CACHE !!! */
197     {
198     unsigned long i, numPreComp;
199     unsigned j;
200     PhotonIdx pIdx;
201     Photon photon;
202     RAY ray;
203     PhotonMap nuPmap;
204 greg 2.1
205 rschregle 2.11 repComplete = numPreComp = finalGather * pmap -> numPhotons;
206 greg 2.1
207 rschregle 2.13 if (verbose) {
208     sprintf(errmsg,
209     "\nPrecomputing irradiance for %ld global photons\n",
210 rschregle 2.11 numPreComp);
211 greg 2.1 eputs(errmsg);
212 rschregle 2.13 #if NIX
213 greg 2.1 fflush(stderr);
214 rschregle 2.13 #endif
215 greg 2.1 }
216    
217 rschregle 2.11 /* Copy photon map for precomputed photons */
218     memcpy(&nuPmap, pmap, sizeof(PhotonMap));
219    
220     /* Zero counters, init new heap and extents */
221     nuPmap.numPhotons = 0;
222     initPhotonHeap(&nuPmap);
223    
224     for (j = 0; j < 3; j++) {
225     nuPmap.minPos [j] = FHUGE;
226     nuPmap.maxPos [j] = -FHUGE;
227 greg 2.1 }
228 rschregle 2.11
229 greg 2.1 /* Record start time, baby */
230     repStartTime = time(NULL);
231 rschregle 2.11 #ifdef SIGCONT
232     signal(SIGCONT, pmapPreCompReport);
233     #endif
234 greg 2.1 repProgress = 0;
235    
236 rschregle 2.11 photonRay(NULL, &ray, PRIMARY, NULL);
237     ray.ro = NULL;
238    
239     for (i = 0; i < numPreComp; i++) {
240     /* Get random photon from stratified distribution in source heap to
241 rschregle 2.13 * avoid duplicates and clustering */
242 rschregle 2.11 pIdx = firstPhoton(pmap) +
243     (unsigned long)((i + pmapRandom(pmap -> randState)) /
244     finalGather);
245     getPhoton(pmap, pIdx, &photon);
246    
247     /* Init dummy photon ray with intersection at photon position */
248     VCOPY(ray.rop, photon.pos);
249     for (j = 0; j < 3; j++)
250     ray.ron [j] = photon.norm [j] / 127.0;
251    
252     /* Get density estimate at photon position */
253     photonDensity(pmap, &ray, ray.rcol);
254    
255     /* Append photon to new heap from ray */
256     newPhoton(&nuPmap, &ray);
257 greg 2.1
258 rschregle 2.11 /* Update progress */
259 greg 2.1 repProgress++;
260    
261     if (photonRepTime > 0 && time(NULL) >= repLastTime + photonRepTime)
262     pmapPreCompReport();
263 rschregle 2.11 #ifdef SIGCONT
264     else signal(SIGCONT, pmapPreCompReport);
265     #endif
266 greg 2.1 }
267    
268 rschregle 2.11 /* Flush heap */
269     flushPhotonHeap(&nuPmap);
270    
271     #ifdef SIGCONT
272     signal(SIGCONT, SIG_DFL);
273     #endif
274    
275     /* Trash original pmap, replace with precomputed one */
276     deletePhotons(pmap);
277     memcpy(pmap, &nuPmap, sizeof(PhotonMap));
278 greg 2.1
279 rschregle 2.13 if (verbose) {
280     eputs("\nRebuilding precomputed photon map\n");
281     #if NIX
282 greg 2.1 fflush(stderr);
283 rschregle 2.13 #endif
284 greg 2.1 }
285 rschregle 2.11
286     /* Rebuild underlying data structure, destroying heap */
287     buildPhotonMap(pmap, NULL, NULL, 1);
288 greg 2.1 }
289    
290    
291    
292 rschregle 2.11 typedef struct {
293     unsigned long numPhotons [NUM_PMAP_TYPES],
294     numEmitted, numComplete;
295     } PhotonCnt;
296    
297    
298    
299     void distribPhotons (PhotonMap **pmaps, unsigned numProc)
300     {
301     EmissionMap emap;
302 rschregle 2.13 char errmsg2 [128], shmFname [PMAP_TMPFNLEN];
303 rschregle 2.11 unsigned t, srcIdx, proc;
304     double totalFlux = 0;
305     int shmFile, stat, pid;
306     PhotonMap *pm;
307     PhotonCnt *photonCnt;
308 greg 2.1
309 rschregle 2.8 for (t = 0; t < NUM_PMAP_TYPES && !pmaps [t]; t++);
310 rschregle 2.11
311 greg 2.1 if (t >= NUM_PMAP_TYPES)
312 rschregle 2.11 error(USER, "no photon maps defined in distribPhotons");
313 greg 2.1
314     if (!nsources)
315 rschregle 2.11 error(USER, "no light sources in distribPhotons");
316 greg 2.1
317     /* ===================================================================
318     * INITIALISATION - Set up emission and scattering funcs
319     * =================================================================== */
320     emap.samples = NULL;
321     emap.maxPartitions = MAXSPART;
322     emap.partitions = (unsigned char*)malloc(emap.maxPartitions >> 1);
323     if (!emap.partitions)
324 rschregle 2.11 error(INTERNAL, "can't allocate source partitions in distribPhotons");
325 greg 2.1
326     /* Initialise all defined photon maps */
327     for (t = 0; t < NUM_PMAP_TYPES; t++)
328 rschregle 2.11 if (pmaps [t]) {
329     initPhotonMap(pmaps [t], t);
330     /* Open photon heapfile */
331     initPhotonHeap(pmaps [t]);
332     /* Per-subprocess target count */
333     pmaps [t] -> distribTarget /= numProc;
334 rschregle 2.13
335     if (!pmaps [t] -> distribTarget)
336     error(INTERNAL, "no photons to distribute in distribPhotons");
337 rschregle 2.11 }
338 greg 2.1
339     initPhotonEmissionFuncs();
340     initPhotonScatterFuncs();
341    
342 rschregle 2.14 /* Get photon ports from modifier list */
343     getPhotonPorts(photonPortList);
344 greg 2.1
345     /* Get photon sensor modifiers */
346     getPhotonSensors(photonSensorList);
347    
348 rschregle 2.13 #if NIX
349 rschregle 2.11 /* Set up shared mem for photon counters (zeroed by ftruncate) */
350 rschregle 2.13 strcpy(shmFname, PMAP_TMPFNAME);
351 rschregle 2.11 shmFile = mkstemp(shmFname);
352    
353 rschregle 2.13 if (shmFile < 0 || ftruncate(shmFile, sizeof(*photonCnt)) < 0)
354     error(SYSTEM, "failed shared mem init in distribPhotons");
355 rschregle 2.11
356     photonCnt = mmap(NULL, sizeof(*photonCnt), PROT_READ | PROT_WRITE,
357     MAP_SHARED, shmFile, 0);
358    
359     if (photonCnt == MAP_FAILED)
360 rschregle 2.13 error(SYSTEM, "failed mapping shared memory in distribPhotons");
361     #else
362     /* Allocate photon counters statically on Windoze */
363     if (!(photonCnt = malloc(sizeof(PhotonCnt))))
364     error(SYSTEM, "failed trivial malloc in distribPhotons");
365     photonCnt -> numEmitted = photonCnt -> numComplete = 0;
366     #endif /* NIX */
367    
368     if (verbose) {
369     sprintf(errmsg, "\nIntegrating flux from %d sources", nsources);
370    
371     if (photonPorts) {
372     sprintf(errmsg2, " via %d ports", numPhotonPorts);
373     strcat(errmsg, errmsg2);
374     }
375    
376     strcat(errmsg, "\n");
377     eputs(errmsg);
378     }
379 greg 2.1
380     /* ===================================================================
381     * FLUX INTEGRATION - Get total photon flux from light sources
382     * =================================================================== */
383 rschregle 2.11 for (srcIdx = 0; srcIdx < nsources; srcIdx++) {
384 greg 2.1 unsigned portCnt = 0;
385     emap.src = source + srcIdx;
386    
387 rschregle 2.11 do { /* Need at least one iteration if no ports! */
388 greg 2.1 emap.port = emap.src -> sflags & SDISTANT ? photonPorts + portCnt
389     : NULL;
390     photonPartition [emap.src -> so -> otype] (&emap);
391    
392 rschregle 2.13 if (verbose) {
393     sprintf(errmsg, "\tIntegrating flux from source %s ",
394 greg 2.1 source [srcIdx].so -> oname);
395 rschregle 2.13
396 greg 2.1 if (emap.port) {
397     sprintf(errmsg2, "via port %s ",
398     photonPorts [portCnt].so -> oname);
399     strcat(errmsg, errmsg2);
400     }
401 rschregle 2.13
402     sprintf(errmsg2, "(%lu partitions)\n", emap.numPartitions);
403 greg 2.1 strcat(errmsg, errmsg2);
404     eputs(errmsg);
405 rschregle 2.13 #if NIX
406 greg 2.1 fflush(stderr);
407 rschregle 2.13 #endif
408 greg 2.1 }
409    
410     for (emap.partitionCnt = 0; emap.partitionCnt < emap.numPartitions;
411     emap.partitionCnt++) {
412     initPhotonEmission(&emap, pdfSamples);
413     totalFlux += colorAvg(emap.partFlux);
414     }
415    
416     portCnt++;
417     } while (portCnt < numPhotonPorts);
418     }
419    
420     if (totalFlux < FTINY)
421     error(USER, "zero flux from light sources");
422 rschregle 2.13
423     /* Record start time for progress reports */
424     repStartTime = time(NULL);
425    
426     if (verbose) {
427     sprintf(errmsg, "\nPhoton distribution @ %d procs\n", numProc);
428     eputs(errmsg);
429     }
430 greg 2.1
431 rschregle 2.11 /* MAIN LOOP */
432     for (proc = 0; proc < numProc; proc++) {
433 rschregle 2.13 #if NIX
434 rschregle 2.11 if (!(pid = fork())) {
435 rschregle 2.13 /* SUBPROCESS ENTERS HERE; open and mmapped files inherited */
436     #else
437     if (1) {
438     /* No subprocess under Windoze */
439     #endif
440     /* Local photon counters for this subprocess */
441 rschregle 2.11 unsigned passCnt = 0, prePassCnt = 0;
442     unsigned long lastNumPhotons [NUM_PMAP_TYPES];
443     unsigned long localNumEmitted = 0; /* Num photons emitted by this
444     subprocess alone */
445 greg 2.1
446 rschregle 2.11 /* Seed RNGs from PID for decorellated photon distribution */
447     pmapSeed(randSeed + proc, partState);
448 rschregle 2.13 pmapSeed(randSeed + (proc + 1) % numProc, emitState);
449     pmapSeed(randSeed + (proc + 2) % numProc, cntState);
450     pmapSeed(randSeed + (proc + 3) % numProc, mediumState);
451     pmapSeed(randSeed + (proc + 4) % numProc, scatterState);
452     pmapSeed(randSeed + (proc + 5) % numProc, rouletteState);
453 rschregle 2.15
454     #ifdef DEBUG_PMAP
455     /* Output child process PID after random delay to prevent corrupted
456     * console output due to race condition */
457     usleep(1e6 * pmapRandom(rouletteState));
458     fprintf(stderr, "Proc %d: PID = %d "
459     "(waiting 10 sec to attach debugger...)\n",
460     proc, getpid());
461     /* Allow time for debugger to attach to child process */
462     sleep(10);
463     #endif
464    
465 greg 2.1 for (t = 0; t < NUM_PMAP_TYPES; t++)
466 rschregle 2.11 lastNumPhotons [t] = 0;
467    
468     /* =============================================================
469     * 2-PASS PHOTON DISTRIBUTION
470     * Pass 1 (pre): emit fraction of target photon count
471     * Pass 2 (main): based on outcome of pass 1, estimate remaining
472     * number of photons to emit to approximate target
473     * count
474     * ============================================================= */
475 greg 2.1 do {
476 rschregle 2.11 double numEmit;
477 greg 2.1
478 rschregle 2.11 if (!passCnt) {
479     /* INIT PASS 1 */
480     /* Skip if no photons contributed after sufficient
481     * iterations; make it clear to user which photon maps are
482     * missing so (s)he can check geometry and materials */
483     if (++prePassCnt > maxPreDistrib) {
484 rschregle 2.13 sprintf(errmsg, "proc %d: too many prepasses", proc);
485 rschregle 2.11
486     for (t = 0; t < NUM_PMAP_TYPES; t++)
487     if (pmaps [t] && !pmaps [t] -> numPhotons) {
488     sprintf(errmsg2, ", no %s photons stored",
489     pmapName [t]);
490     strcat(errmsg, errmsg2);
491     }
492    
493     error(USER, errmsg);
494     break;
495 greg 2.1 }
496 rschregle 2.11
497     /* Num to emit is fraction of minimum target count */
498     numEmit = FHUGE;
499 greg 2.1
500 rschregle 2.11 for (t = 0; t < NUM_PMAP_TYPES; t++)
501     if (pmaps [t])
502     numEmit = min(pmaps [t] -> distribTarget, numEmit);
503    
504     numEmit *= preDistrib;
505 greg 2.1 }
506 rschregle 2.11 else {
507     /* INIT PASS 2 */
508     /* Based on the outcome of the predistribution we can now
509     * estimate how many more photons we have to emit for each
510     * photon map to meet its respective target count. This
511     * value is clamped to 0 in case the target has already been
512     * exceeded in the pass 1. */
513     double maxDistribRatio = 0;
514    
515     /* Set the distribution ratio for each map; this indicates
516     * how many photons of each respective type are stored per
517     * emitted photon, and is used as probability for storing a
518     * photon by newPhoton(). Since this biases the photon
519     * density, newPhoton() promotes the flux of stored photons
520     * to compensate. */
521     for (t = 0; t < NUM_PMAP_TYPES; t++)
522     if ((pm = pmaps [t])) {
523     pm -> distribRatio = (double)pm -> distribTarget /
524     pm -> numPhotons - 1;
525    
526     /* Check if photon map "overflowed", i.e. exceeded its
527     * target count in the prepass; correcting the photon
528     * flux via the distribution ratio is no longer
529     * possible, as no more photons of this type will be
530     * stored, so notify the user rather than deliver
531     * incorrect results. In future we should handle this
532     * more intelligently by using the photonFlux in each
533     * photon map to individually correct the flux after
534     * distribution. */
535     if (pm -> distribRatio <= FTINY) {
536     sprintf(errmsg, "%s photon map overflow in "
537     "prepass, reduce -apD", pmapName [t]);
538     error(INTERNAL, errmsg);
539     }
540    
541     maxDistribRatio = max(pm -> distribRatio,
542     maxDistribRatio);
543     }
544 greg 2.1
545 rschregle 2.11 /* Normalise distribution ratios and calculate number of
546     * photons to emit in main pass */
547     for (t = 0; t < NUM_PMAP_TYPES; t++)
548     if ((pm = pmaps [t]))
549     pm -> distribRatio /= maxDistribRatio;
550    
551     if ((numEmit = localNumEmitted * maxDistribRatio) < FTINY)
552     /* No photons left to distribute in main pass */
553     break;
554     }
555    
556 rschregle 2.13 /* Update shared completion counter for progreport by parent */
557 rschregle 2.11 photonCnt -> numComplete += numEmit;
558    
559     /* PHOTON DISTRIBUTION LOOP */
560     for (srcIdx = 0; srcIdx < nsources; srcIdx++) {
561     unsigned portCnt = 0;
562     emap.src = source + srcIdx;
563    
564     do { /* Need at least one iteration if no ports! */
565     emap.port = emap.src -> sflags & SDISTANT
566     ? photonPorts + portCnt : NULL;
567     photonPartition [emap.src -> so -> otype] (&emap);
568    
569 rschregle 2.13 if (verbose && !proc) {
570     /* Output from subproc 0 only to avoid race condition
571     * on console I/O */
572 rschregle 2.11 if (!passCnt)
573 rschregle 2.13 sprintf(errmsg, "\tPREPASS %d on source %s ",
574 rschregle 2.11 prePassCnt, source [srcIdx].so -> oname);
575     else
576 rschregle 2.13 sprintf(errmsg, "\tMAIN PASS on source %s ",
577 rschregle 2.11 source [srcIdx].so -> oname);
578 rschregle 2.13
579 rschregle 2.11 if (emap.port) {
580     sprintf(errmsg2, "via port %s ",
581     photonPorts [portCnt].so -> oname);
582     strcat(errmsg, errmsg2);
583     }
584 rschregle 2.13
585 rschregle 2.11 sprintf(errmsg2, "(%lu partitions)\n",
586     emap.numPartitions);
587     strcat(errmsg, errmsg2);
588     eputs(errmsg);
589 rschregle 2.13 #if NIX
590 rschregle 2.11 fflush(stderr);
591 rschregle 2.13 #endif
592 rschregle 2.11 }
593 greg 2.1
594 rschregle 2.11 for (emap.partitionCnt = 0; emap.partitionCnt < emap.numPartitions;
595     emap.partitionCnt++) {
596     double partNumEmit;
597     unsigned long partEmitCnt;
598    
599     /* Get photon origin within current source partishunn
600     * and build emission map */
601     photonOrigin [emap.src -> so -> otype] (&emap);
602     initPhotonEmission(&emap, pdfSamples);
603    
604     /* Number of photons to emit from ziss partishunn --
605     * proportional to flux; photon ray weight and scalar
606 rschregle 2.13 * flux are uniform (latter only varying in RGB). */
607 rschregle 2.11 partNumEmit = numEmit * colorAvg(emap.partFlux) /
608     totalFlux;
609     partEmitCnt = (unsigned long)partNumEmit;
610    
611     /* Probabilistically account for fractional photons */
612     if (pmapRandom(cntState) < partNumEmit - partEmitCnt)
613     partEmitCnt++;
614    
615     /* Update local and shared (global) emission counter */
616     photonCnt -> numEmitted += partEmitCnt;
617     localNumEmitted += partEmitCnt;
618    
619     /* Integer counter avoids FP rounding errors during
620     * iteration */
621     while (partEmitCnt--) {
622     RAY photonRay;
623    
624     /* Emit photon based on PDF and trace through scene
625     * until absorbed/leaked */
626     emitPhoton(&emap, &photonRay);
627 rschregle 2.13 #if 1
628     if (emap.port)
629     /* !!! PHOTON PORT REJECTION SAMPLING HACK: set
630     * !!! photon port as fake hit object for
631     * !!! primary ray to check for intersection in
632     * !!! tracePhoton() */
633     photonRay.ro = emap.port -> so;
634     #endif
635 rschregle 2.11 tracePhoton(&photonRay);
636     }
637 rschregle 2.13
638 rschregle 2.11 /* Update shared global photon count for each pmap */
639     for (t = 0; t < NUM_PMAP_TYPES; t++)
640     if (pmaps [t]) {
641     photonCnt -> numPhotons [t] +=
642     pmaps [t] -> numPhotons - lastNumPhotons [t];
643     lastNumPhotons [t] = pmaps [t] -> numPhotons;
644     }
645 rschregle 2.13 #if !NIX
646     /* Synchronous progress report on Windoze */
647     if (!proc && photonRepTime > 0 &&
648     time(NULL) >= repLastTime + photonRepTime) {
649     repEmitted = repProgress = photonCnt -> numEmitted;
650     repComplete = photonCnt -> numComplete;
651     pmapDistribReport();
652     }
653     #endif
654 rschregle 2.11 }
655 greg 2.1
656 rschregle 2.11 portCnt++;
657     } while (portCnt < numPhotonPorts);
658     }
659    
660     for (t = 0; t < NUM_PMAP_TYPES; t++)
661     if (pmaps [t] && !pmaps [t] -> numPhotons) {
662     /* Double preDistrib in case a photon map is empty and
663     * redo pass 1 --> possibility of infinite loop for
664     * pathological scenes (e.g. absorbing materials) */
665     preDistrib *= 2;
666     break;
667 greg 2.1 }
668 rschregle 2.11
669 rschregle 2.13 if (t >= NUM_PMAP_TYPES)
670 rschregle 2.11 /* No empty photon maps found; now do pass 2 */
671     passCnt++;
672     } while (passCnt < 2);
673    
674 rschregle 2.13 /* Flush heap buffa for every pmap one final time;
675     * avoids potential data corruption! */
676 rschregle 2.11 for (t = 0; t < NUM_PMAP_TYPES; t++)
677     if (pmaps [t]) {
678     flushPhotonHeap(pmaps [t]);
679 rschregle 2.13 /* Heap file closed automatically on exit
680     fclose(pmaps [t] -> heap); */
681 rschregle 2.11 #ifdef DEBUG_PMAP
682 rschregle 2.13 sprintf(errmsg, "Proc %d: total %ld photons\n", proc,
683 rschregle 2.11 pmaps [t] -> numPhotons);
684     eputs(errmsg);
685     #endif
686     }
687 rschregle 2.13 #if NIX
688     /* Terminate subprocess */
689 rschregle 2.11 exit(0);
690 rschregle 2.13 #endif
691 greg 2.1 }
692 rschregle 2.11 else if (pid < 0)
693     error(SYSTEM, "failed to fork subprocess in distribPhotons");
694     }
695    
696 rschregle 2.13 #if NIX
697 rschregle 2.11 /* PARENT PROCESS CONTINUES HERE */
698     #ifdef SIGCONT
699 rschregle 2.13 /* Enable progress report signal handler */
700 rschregle 2.11 signal(SIGCONT, pmapDistribReport);
701 rschregle 2.13 #endif
702     /* Wait for subprocesses complete while reporting progress */
703 rschregle 2.11 proc = numProc;
704     while (proc) {
705     while (waitpid(-1, &stat, WNOHANG) > 0) {
706     /* Subprocess exited; check status */
707     if (!WIFEXITED(stat) || WEXITSTATUS(stat))
708     error(USER, "failed photon distribution");
709    
710     --proc;
711     }
712    
713     /* Nod off for a bit and update progress */
714     sleep(1);
715 rschregle 2.13
716     /* Asynchronous progress report from shared subprocess counters */
717 rschregle 2.11 repEmitted = repProgress = photonCnt -> numEmitted;
718 rschregle 2.13 repComplete = photonCnt -> numComplete;
719 rschregle 2.11
720 rschregle 2.13 repProgress = repComplete = 0;
721 greg 2.1 for (t = 0; t < NUM_PMAP_TYPES; t++)
722 rschregle 2.11 if ((pm = pmaps [t])) {
723     /* Get global photon count from shmem updated by subprocs */
724 rschregle 2.13 repProgress += pm -> numPhotons = photonCnt -> numPhotons [t];
725     repComplete += pm -> distribTarget;
726 greg 2.1 }
727 rschregle 2.13 repComplete *= numProc;
728 rschregle 2.11
729     if (photonRepTime > 0 && time(NULL) >= repLastTime + photonRepTime)
730     pmapDistribReport();
731     #ifdef SIGCONT
732     else signal(SIGCONT, pmapDistribReport);
733     #endif
734     }
735 rschregle 2.13 #endif /* NIX */
736 greg 2.1
737     /* ===================================================================
738 rschregle 2.11 * POST-DISTRIBUTION - Set photon flux and build data struct for photon
739     * storage, etc.
740 greg 2.1 * =================================================================== */
741 rschregle 2.11 #ifdef SIGCONT
742 rschregle 2.13 /* Reset signal handler */
743 rschregle 2.11 signal(SIGCONT, SIG_DFL);
744     #endif
745 greg 2.1 free(emap.samples);
746    
747 rschregle 2.13 /* Set photon flux */
748 rschregle 2.11 totalFlux /= photonCnt -> numEmitted;
749 rschregle 2.13 #if NIX
750 rschregle 2.11 /* Photon counters no longer needed, unmap shared memory */
751     munmap(photonCnt, sizeof(*photonCnt));
752     close(shmFile);
753 rschregle 2.13 unlink(shmFname);
754 rschregle 2.11 #else
755 rschregle 2.13 free(photonCnt);
756 rschregle 2.11 #endif
757 rschregle 2.13 if (verbose)
758     eputs("\n");
759    
760 greg 2.1 for (t = 0; t < NUM_PMAP_TYPES; t++)
761 rschregle 2.8 if (pmaps [t]) {
762 rschregle 2.13 if (verbose) {
763     sprintf(errmsg, "Building %s photon map\n", pmapName [t]);
764 greg 2.1 eputs(errmsg);
765 rschregle 2.13 #if NIX
766 greg 2.1 fflush(stderr);
767 rschregle 2.13 #endif
768 greg 2.1 }
769 rschregle 2.11
770     /* Build underlying data structure; heap is destroyed */
771     buildPhotonMap(pmaps [t], &totalFlux, NULL, numProc);
772 greg 2.1 }
773 rschregle 2.13
774 greg 2.1 /* Precompute photon irradiance if necessary */
775 rschregle 2.13 if (preCompPmap) {
776     if (verbose)
777     eputs("\n");
778 greg 2.1 preComputeGlobal(preCompPmap);
779 rschregle 2.13 }
780    
781     if (verbose)
782     eputs("\n");
783 greg 2.1 }