ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/rt/pmap.c
Revision: 2.16
Committed: Wed Nov 21 19:33:30 2018 UTC (6 years, 5 months ago) by rschregle
Content type: text/plain
Branch: MAIN
Changes since 2.15: +24 -14 lines
Log Message:
Fixed potential infinite loop in photonParticipate() if photon leaves scene

File Contents

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