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 (5 years, 11 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

# Content
1 #ifndef lint
2 static const char RCSid[] = "$Id: pmap.c,v 2.14 2018/03/20 19:55:33 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.14 2018/03/20 19:55:33 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 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 #if 0
73 if (albedo > FTINY && ray -> rlvl > 0)
74 #else
75 /* Store volume photons unconditionally in mist to also account for
76 direct inscattering from sources */
77 if (albedo > FTINY)
78 #endif
79 /* Add to volume photon map */
80 newPhoton(volumePmap, ray);
81
82 /* Absorbed? */
83 if (pmapRandom(rouletteState) > albedo)
84 return 0;
85
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 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
147 if (ray -> rlvl > photonMaxBounce) {
148 #ifdef PMAP_RUNAWAY_WARN
149 error(WARNING, "runaway photon!");
150 #endif
151 return;
152 }
153
154 if (colorAvg(ray -> cext) > FTINY && !photonParticipate(ray))
155 return;
156
157 if (localhit(ray, &thescene)) {
158 mod = ray -> ro -> omod;
159
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 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 /* 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
205 repComplete = numPreComp = finalGather * pmap -> numPhotons;
206
207 if (verbose) {
208 sprintf(errmsg,
209 "\nPrecomputing irradiance for %ld global photons\n",
210 numPreComp);
211 eputs(errmsg);
212 #if NIX
213 fflush(stderr);
214 #endif
215 }
216
217 /* 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 }
228
229 /* Record start time, baby */
230 repStartTime = time(NULL);
231 #ifdef SIGCONT
232 signal(SIGCONT, pmapPreCompReport);
233 #endif
234 repProgress = 0;
235
236 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 * avoid duplicates and clustering */
242 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
258 /* Update progress */
259 repProgress++;
260
261 if (photonRepTime > 0 && time(NULL) >= repLastTime + photonRepTime)
262 pmapPreCompReport();
263 #ifdef SIGCONT
264 else signal(SIGCONT, pmapPreCompReport);
265 #endif
266 }
267
268 /* 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
279 if (verbose) {
280 eputs("\nRebuilding precomputed photon map\n");
281 #if NIX
282 fflush(stderr);
283 #endif
284 }
285
286 /* Rebuild underlying data structure, destroying heap */
287 buildPhotonMap(pmap, NULL, NULL, 1);
288 }
289
290
291
292 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 char errmsg2 [128], shmFname [PMAP_TMPFNLEN];
303 unsigned t, srcIdx, proc;
304 double totalFlux = 0;
305 int shmFile, stat, pid;
306 PhotonMap *pm;
307 PhotonCnt *photonCnt;
308
309 for (t = 0; t < NUM_PMAP_TYPES && !pmaps [t]; t++);
310
311 if (t >= NUM_PMAP_TYPES)
312 error(USER, "no photon maps defined in distribPhotons");
313
314 if (!nsources)
315 error(USER, "no light sources in distribPhotons");
316
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 error(INTERNAL, "can't allocate source partitions in distribPhotons");
325
326 /* Initialise all defined photon maps */
327 for (t = 0; t < NUM_PMAP_TYPES; t++)
328 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
335 if (!pmaps [t] -> distribTarget)
336 error(INTERNAL, "no photons to distribute in distribPhotons");
337 }
338
339 initPhotonEmissionFuncs();
340 initPhotonScatterFuncs();
341
342 /* Get photon ports from modifier list */
343 getPhotonPorts(photonPortList);
344
345 /* Get photon sensor modifiers */
346 getPhotonSensors(photonSensorList);
347
348 #if NIX
349 /* Set up shared mem for photon counters (zeroed by ftruncate) */
350 strcpy(shmFname, PMAP_TMPFNAME);
351 shmFile = mkstemp(shmFname);
352
353 if (shmFile < 0 || ftruncate(shmFile, sizeof(*photonCnt)) < 0)
354 error(SYSTEM, "failed shared mem init in distribPhotons");
355
356 photonCnt = mmap(NULL, sizeof(*photonCnt), PROT_READ | PROT_WRITE,
357 MAP_SHARED, shmFile, 0);
358
359 if (photonCnt == MAP_FAILED)
360 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
380 /* ===================================================================
381 * FLUX INTEGRATION - Get total photon flux from light sources
382 * =================================================================== */
383 for (srcIdx = 0; srcIdx < nsources; srcIdx++) {
384 unsigned portCnt = 0;
385 emap.src = source + srcIdx;
386
387 do { /* Need at least one iteration if no ports! */
388 emap.port = emap.src -> sflags & SDISTANT ? photonPorts + portCnt
389 : NULL;
390 photonPartition [emap.src -> so -> otype] (&emap);
391
392 if (verbose) {
393 sprintf(errmsg, "\tIntegrating flux from source %s ",
394 source [srcIdx].so -> oname);
395
396 if (emap.port) {
397 sprintf(errmsg2, "via port %s ",
398 photonPorts [portCnt].so -> oname);
399 strcat(errmsg, errmsg2);
400 }
401
402 sprintf(errmsg2, "(%lu partitions)\n", emap.numPartitions);
403 strcat(errmsg, errmsg2);
404 eputs(errmsg);
405 #if NIX
406 fflush(stderr);
407 #endif
408 }
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
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
431 /* MAIN LOOP */
432 for (proc = 0; proc < numProc; proc++) {
433 #if NIX
434 if (!(pid = fork())) {
435 /* 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 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
446 /* Seed RNGs from PID for decorellated photon distribution */
447 pmapSeed(randSeed + proc, partState);
448 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
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 for (t = 0; t < NUM_PMAP_TYPES; t++)
466 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 do {
476 double numEmit;
477
478 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 sprintf(errmsg, "proc %d: too many prepasses", proc);
485
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 }
496
497 /* Num to emit is fraction of minimum target count */
498 numEmit = FHUGE;
499
500 for (t = 0; t < NUM_PMAP_TYPES; t++)
501 if (pmaps [t])
502 numEmit = min(pmaps [t] -> distribTarget, numEmit);
503
504 numEmit *= preDistrib;
505 }
506 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
545 /* 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 /* Update shared completion counter for progreport by parent */
557 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 if (verbose && !proc) {
570 /* Output from subproc 0 only to avoid race condition
571 * on console I/O */
572 if (!passCnt)
573 sprintf(errmsg, "\tPREPASS %d on source %s ",
574 prePassCnt, source [srcIdx].so -> oname);
575 else
576 sprintf(errmsg, "\tMAIN PASS on source %s ",
577 source [srcIdx].so -> oname);
578
579 if (emap.port) {
580 sprintf(errmsg2, "via port %s ",
581 photonPorts [portCnt].so -> oname);
582 strcat(errmsg, errmsg2);
583 }
584
585 sprintf(errmsg2, "(%lu partitions)\n",
586 emap.numPartitions);
587 strcat(errmsg, errmsg2);
588 eputs(errmsg);
589 #if NIX
590 fflush(stderr);
591 #endif
592 }
593
594 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 * flux are uniform (latter only varying in RGB). */
607 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 #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 tracePhoton(&photonRay);
636 }
637
638 /* 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 #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 }
655
656 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 }
668
669 if (t >= NUM_PMAP_TYPES)
670 /* No empty photon maps found; now do pass 2 */
671 passCnt++;
672 } while (passCnt < 2);
673
674 /* Flush heap buffa for every pmap one final time;
675 * avoids potential data corruption! */
676 for (t = 0; t < NUM_PMAP_TYPES; t++)
677 if (pmaps [t]) {
678 flushPhotonHeap(pmaps [t]);
679 /* Heap file closed automatically on exit
680 fclose(pmaps [t] -> heap); */
681 #ifdef DEBUG_PMAP
682 sprintf(errmsg, "Proc %d: total %ld photons\n", proc,
683 pmaps [t] -> numPhotons);
684 eputs(errmsg);
685 #endif
686 }
687 #if NIX
688 /* Terminate subprocess */
689 exit(0);
690 #endif
691 }
692 else if (pid < 0)
693 error(SYSTEM, "failed to fork subprocess in distribPhotons");
694 }
695
696 #if NIX
697 /* PARENT PROCESS CONTINUES HERE */
698 #ifdef SIGCONT
699 /* Enable progress report signal handler */
700 signal(SIGCONT, pmapDistribReport);
701 #endif
702 /* Wait for subprocesses complete while reporting progress */
703 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
716 /* Asynchronous progress report from shared subprocess counters */
717 repEmitted = repProgress = photonCnt -> numEmitted;
718 repComplete = photonCnt -> numComplete;
719
720 repProgress = repComplete = 0;
721 for (t = 0; t < NUM_PMAP_TYPES; t++)
722 if ((pm = pmaps [t])) {
723 /* Get global photon count from shmem updated by subprocs */
724 repProgress += pm -> numPhotons = photonCnt -> numPhotons [t];
725 repComplete += pm -> distribTarget;
726 }
727 repComplete *= numProc;
728
729 if (photonRepTime > 0 && time(NULL) >= repLastTime + photonRepTime)
730 pmapDistribReport();
731 #ifdef SIGCONT
732 else signal(SIGCONT, pmapDistribReport);
733 #endif
734 }
735 #endif /* NIX */
736
737 /* ===================================================================
738 * POST-DISTRIBUTION - Set photon flux and build data struct for photon
739 * storage, etc.
740 * =================================================================== */
741 #ifdef SIGCONT
742 /* Reset signal handler */
743 signal(SIGCONT, SIG_DFL);
744 #endif
745 free(emap.samples);
746
747 /* Set photon flux */
748 totalFlux /= photonCnt -> numEmitted;
749 #if NIX
750 /* Photon counters no longer needed, unmap shared memory */
751 munmap(photonCnt, sizeof(*photonCnt));
752 close(shmFile);
753 unlink(shmFname);
754 #else
755 free(photonCnt);
756 #endif
757 if (verbose)
758 eputs("\n");
759
760 for (t = 0; t < NUM_PMAP_TYPES; t++)
761 if (pmaps [t]) {
762 if (verbose) {
763 sprintf(errmsg, "Building %s photon map\n", pmapName [t]);
764 eputs(errmsg);
765 #if NIX
766 fflush(stderr);
767 #endif
768 }
769
770 /* Build underlying data structure; heap is destroyed */
771 buildPhotonMap(pmaps [t], &totalFlux, NULL, numProc);
772 }
773
774 /* Precompute photon irradiance if necessary */
775 if (preCompPmap) {
776 if (verbose)
777 eputs("\n");
778 preComputeGlobal(preCompPmap);
779 }
780
781 if (verbose)
782 eputs("\n");
783 }