ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/rt/pmapcontrib.c
Revision: 2.16
Committed: Fri Mar 16 21:00:09 2018 UTC (7 years, 1 month ago) by rschregle
Content type: text/plain
Branch: MAIN
Changes since 2.15: +79 -18 lines
Log Message:
Fixed serious endless predistrib loop bug in distribPhotonContrib(), added
on-demand debugging info via SIGUSR1 from child processes.

File Contents

# User Rev Content
1 greg 2.9 #ifndef lint
2 rschregle 2.16 static const char RCSid[] = "$Id: pmapcontrib.c,v 2.15 2018/02/08 19:55:02 rschregle Exp $";
3 greg 2.9 #endif
4 rschregle 2.12
5 greg 2.1 /*
6 rschregle 2.12 ======================================================================
7 rschregle 2.14 Photon map for light source contributions
8 greg 2.1
9     Roland Schregle (roland.schregle@{hslu.ch, gmail.com})
10 rschregle 2.4 (c) Lucerne University of Applied Sciences and Arts,
11 rschregle 2.12 supported by the Swiss National Science Foundation (SNSF, #147053)
12     ======================================================================
13 greg 2.1
14 rschregle 2.16 $Id: pmapcontrib.c,v 2.15 2018/02/08 19:55:02 rschregle Exp $
15 greg 2.1 */
16    
17    
18     #include "pmapcontrib.h"
19     #include "pmapmat.h"
20     #include "pmapsrc.h"
21     #include "pmaprand.h"
22     #include "pmapio.h"
23     #include "pmapdiag.h"
24     #include "rcontrib.h"
25     #include "otypes.h"
26 rschregle 2.14 #if NIX
27     #include <sys/mman.h>
28 rschregle 2.15 #include <sys/wait.h>
29 rschregle 2.14 #endif
30 greg 2.1
31    
32 rschregle 2.12 static PhotonPrimaryIdx newPhotonPrimary (PhotonMap *pmap,
33     const RAY *primRay,
34     FILE *primHeap)
35     /* Add primary ray for emitted photon and save light source index, origin on
36     * source, and emitted direction; used by contrib photons. The current
37     * primary is stored in pmap -> lastPrimary. If the previous primary
38     * contributed photons (has srcIdx >= 0), it's appended to primHeap. If
39     * primRay == NULL, the current primary is still flushed, but no new primary
40     * is set. Returns updated primary counter pmap -> numPrimary. */
41 greg 2.1 {
42 rschregle 2.12 if (!pmap || !primHeap)
43     return 0;
44 greg 2.1
45 rschregle 2.12 /* Check if last primary ray has spawned photons (srcIdx >= 0, see
46 rschregle 2.14 * newPhoton()), in which case we save it to the primary heap file
47     * before clobbering it */
48 rschregle 2.12 if (pmap -> lastPrimary.srcIdx >= 0) {
49     if (!fwrite(&pmap -> lastPrimary, sizeof(PhotonPrimary), 1, primHeap))
50     error(SYSTEM, "failed writing photon primary in newPhotonPrimary");
51    
52     pmap -> numPrimary++;
53     if (pmap -> numPrimary > PMAP_MAXPRIMARY)
54     error(INTERNAL, "photon primary overflow in newPhotonPrimary");
55     }
56    
57     /* Mark unused with negative source index until path spawns a photon (see
58     * newPhoton()) */
59     pmap -> lastPrimary.srcIdx = -1;
60    
61     if (primRay) {
62     FVECT dvec;
63 rschregle 2.15
64     #ifdef PMAP_PRIMARYDIR
65 rschregle 2.12 /* Reverse incident direction to point to light source */
66     dvec [0] = -primRay -> rdir [0];
67     dvec [1] = -primRay -> rdir [1];
68     dvec [2] = -primRay -> rdir [2];
69     pmap -> lastPrimary.dir = encodedir(dvec);
70 rschregle 2.15 #endif
71 rschregle 2.12 #ifdef PMAP_PRIMARYPOS
72     VCOPY(pmap -> lastPrimary.pos, primRay -> rop);
73     #endif
74     }
75    
76     return pmap -> numPrimary;
77     }
78    
79 greg 2.1
80    
81 rschregle 2.14 #ifdef DEBUG_PMAP
82 rschregle 2.12 static int checkPrimaryHeap (FILE *file)
83     /* Check heap for ordered primaries */
84     {
85     Photon p, lastp;
86     int i, dup;
87    
88     rewind(file);
89     memset(&lastp, 0, sizeof(lastp));
90 greg 2.1
91 rschregle 2.12 while (fread(&p, sizeof(p), 1, file)) {
92     dup = 1;
93    
94     for (i = 0; i <= 2; i++) {
95     if (p.pos [i] < thescene.cuorg [i] ||
96     p.pos [i] > thescene.cuorg [i] + thescene.cusize) {
97    
98     sprintf(errmsg, "corrupt photon in heap at [%f, %f, %f]\n",
99     p.pos [0], p.pos [1], p.pos [2]);
100     error(WARNING, errmsg);
101     }
102    
103     dup &= p.pos [i] == lastp.pos [i];
104     }
105    
106     if (dup) {
107     sprintf(errmsg,
108     "consecutive duplicate photon in heap at [%f, %f, %f]\n",
109     p.pos [0], p.pos [1], p.pos [2]);
110 greg 2.1 error(WARNING, errmsg);
111 rschregle 2.12 }
112 greg 2.1 }
113    
114 rschregle 2.12 return 0;
115     }
116     #endif
117    
118    
119    
120     static PhotonPrimaryIdx buildPrimaries (PhotonMap *pmap, FILE **primaryHeap,
121 rschregle 2.14 char **primaryHeapFname,
122 rschregle 2.12 PhotonPrimaryIdx *primaryOfs,
123     unsigned numHeaps)
124     /* Consolidate per-subprocess photon primary heaps into the primary array
125     * pmap -> primaries. Returns offset for primary index linearisation in
126     * numPrimary. The heap files in primaryHeap are closed on return. */
127     {
128     PhotonPrimaryIdx heapLen;
129     unsigned heap;
130 greg 2.1
131 rschregle 2.12 if (!pmap || !primaryHeap || !primaryOfs || !numHeaps)
132     return 0;
133 greg 2.1
134 rschregle 2.12 pmap -> numPrimary = 0;
135    
136     for (heap = 0; heap < numHeaps; heap++) {
137     primaryOfs [heap] = pmap -> numPrimary;
138    
139 rschregle 2.14 if (fseek(primaryHeap [heap], 0, SEEK_END) < 0)
140 rschregle 2.12 error(SYSTEM, "failed photon primary seek in buildPrimaries");
141     pmap -> numPrimary += heapLen = ftell(primaryHeap [heap]) /
142     sizeof(PhotonPrimary);
143    
144     pmap -> primaries = realloc(pmap -> primaries,
145     pmap -> numPrimary *
146     sizeof(PhotonPrimary));
147     if (!pmap -> primaries)
148     error(SYSTEM, "failed photon primary alloc in buildPrimaries");
149    
150     rewind(primaryHeap [heap]);
151     if (fread(pmap -> primaries + primaryOfs [heap], sizeof(PhotonPrimary),
152     heapLen, primaryHeap [heap]) != heapLen)
153     error(SYSTEM, "failed reading photon primaries in buildPrimaries");
154 greg 2.1
155 rschregle 2.14 fclose(primaryHeap [heap]);
156     unlink(primaryHeapFname [heap]);
157 rschregle 2.12 }
158    
159     return pmap -> numPrimary;
160     }
161 greg 2.6
162    
163 greg 2.7
164 rschregle 2.12 /* Defs for photon emission counter array passed by sub-processes to parent
165     * via shared memory */
166     typedef unsigned long PhotonContribCnt;
167 greg 2.7
168 rschregle 2.12 /* Indices for photon emission counter array: num photons stored and num
169     * emitted per source */
170     #define PHOTONCNT_NUMPHOT 0
171     #define PHOTONCNT_NUMEMIT(n) (1 + n)
172 greg 2.1
173    
174    
175 rschregle 2.16
176    
177    
178 rschregle 2.12 void distribPhotonContrib (PhotonMap* pm, unsigned numProc)
179 greg 2.1 {
180 rschregle 2.12 EmissionMap emap;
181 rschregle 2.14 char errmsg2 [128], shmFname [PMAP_TMPFNLEN];
182 rschregle 2.12 unsigned srcIdx, proc;
183     int shmFile, stat, pid;
184     double *srcFlux, /* Emitted flux per light source */
185     srcDistribTarget; /* Target photon count per source */
186     PhotonContribCnt *photonCnt; /* Photon emission counter array */
187 rschregle 2.14 unsigned photonCntSize = sizeof(PhotonContribCnt) *
188 rschregle 2.12 PHOTONCNT_NUMEMIT(nsources);
189 rschregle 2.14 FILE **primaryHeap = NULL;
190     char **primaryHeapFname = NULL;
191     PhotonPrimaryIdx *primaryOfs = NULL;
192 rschregle 2.12
193 greg 2.1 if (!pm)
194 rschregle 2.12 error(USER, "no photon map defined in distribPhotonContrib");
195 greg 2.1
196     if (!nsources)
197 rschregle 2.12 error(USER, "no light sources in distribPhotonContrib");
198    
199     if (nsources > MAXMODLIST)
200     error(USER, "too many light sources in distribPhotonContrib");
201    
202 greg 2.1 /* Allocate photon flux per light source; this differs for every
203     * source as all sources contribute the same number of distributed
204     * photons (srcDistribTarget), hence the number of photons emitted per
205     * source does not correlate with its emitted flux. The resulting flux
206     * per photon is therefore adjusted individually for each source. */
207     if (!(srcFlux = calloc(nsources, sizeof(double))))
208 rschregle 2.12 error(SYSTEM, "can't allocate source flux in distribPhotonContrib");
209 greg 2.1
210 rschregle 2.12 /* ===================================================================
211     * INITIALISATION - Set up emission and scattering funcs
212     * =================================================================== */
213 greg 2.1 emap.samples = NULL;
214     emap.src = NULL;
215     emap.maxPartitions = MAXSPART;
216     emap.partitions = (unsigned char*)malloc(emap.maxPartitions >> 1);
217     if (!emap.partitions)
218 rschregle 2.12 error(USER, "can't allocate source partitions in distribPhotonContrib");
219 greg 2.1
220 rschregle 2.12 /* Initialise contrib photon map */
221 greg 2.1 initPhotonMap(pm, PMAP_TYPE_CONTRIB);
222 rschregle 2.12 initPhotonHeap(pm);
223 greg 2.1 initPhotonEmissionFuncs();
224     initPhotonScatterFuncs();
225    
226 rschregle 2.12 /* Per-subprocess / per-source target counts */
227     pm -> distribTarget /= numProc;
228 rschregle 2.14 srcDistribTarget = nsources ? (double)pm -> distribTarget / nsources : 0;
229    
230     if (!pm -> distribTarget)
231     error(INTERNAL, "no photons to distribute in distribPhotonContrib");
232 rschregle 2.12
233 greg 2.1 /* Get photon ports if specified */
234     if (ambincl == 1)
235     getPhotonPorts();
236    
237     /* Get photon sensor modifiers */
238     getPhotonSensors(photonSensorList);
239 rschregle 2.14
240     #if NIX
241 rschregle 2.12 /* Set up shared mem for photon counters (zeroed by ftruncate) */
242 rschregle 2.14 strcpy(shmFname, PMAP_TMPFNAME);
243 rschregle 2.12 shmFile = mkstemp(shmFname);
244    
245     if (shmFile < 0 || ftruncate(shmFile, photonCntSize) < 0)
246     error(SYSTEM, "failed shared mem init in distribPhotonContrib");
247 greg 2.1
248 rschregle 2.12 photonCnt = mmap(NULL, photonCntSize, PROT_READ | PROT_WRITE,
249     MAP_SHARED, shmFile, 0);
250    
251     if (photonCnt == MAP_FAILED)
252     error(SYSTEM, "failed shared mem mapping in distribPhotonContrib");
253 rschregle 2.14 #else
254     /* Allocate photon counters statically on Windoze */
255     if (!(photonCnt = malloc(photonCntSize)))
256     error(SYSTEM, "failed trivial malloc in distribPhotonContrib");
257    
258     for (srcIdx = 0; srcIdx < PHOTONCNT_NUMEMIT(nsources); srcIdx++)
259     photonCnt [srcIdx] = 0;
260     #endif /* NIX */
261    
262     if (verbose) {
263     sprintf(errmsg, "\nIntegrating flux from %d sources", nsources);
264    
265     if (photonPorts) {
266     sprintf(errmsg2, " via %d ports", numPhotonPorts);
267     strcat(errmsg, errmsg2);
268     }
269    
270     strcat(errmsg, "\n");
271     eputs(errmsg);
272     }
273 rschregle 2.12
274     /* =============================================================
275 rschregle 2.14 * FLUX INTEGRATION - Get total flux emitted from sources/ports
276 rschregle 2.12 * ============================================================= */
277 greg 2.1 for (srcIdx = 0; srcIdx < nsources; srcIdx++) {
278 rschregle 2.14 unsigned portCnt = 0;
279 rschregle 2.12 srcFlux [srcIdx] = 0;
280 greg 2.1 emap.src = source + srcIdx;
281    
282 rschregle 2.12 do { /* Need at least one iteration if no ports! */
283     emap.port = emap.src -> sflags & SDISTANT ? photonPorts + portCnt
284     : NULL;
285 greg 2.1 photonPartition [emap.src -> so -> otype] (&emap);
286 rschregle 2.14
287     if (verbose) {
288     sprintf(errmsg, "\tIntegrating flux from source %s ",
289     source [srcIdx].so -> oname);
290    
291 greg 2.1 if (emap.port) {
292     sprintf(errmsg2, "via port %s ",
293     photonPorts [portCnt].so -> oname);
294     strcat(errmsg, errmsg2);
295     }
296 rschregle 2.14
297     sprintf(errmsg2, "(%lu partitions)\n", emap.numPartitions);
298 greg 2.1 strcat(errmsg, errmsg2);
299     eputs(errmsg);
300 rschregle 2.14 #if NIX
301 greg 2.1 fflush(stderr);
302 rschregle 2.14 #endif
303     }
304 greg 2.1
305 rschregle 2.12 for (emap.partitionCnt = 0; emap.partitionCnt < emap.numPartitions;
306 greg 2.1 emap.partitionCnt++) {
307     initPhotonEmission(&emap, pdfSamples);
308     srcFlux [srcIdx] += colorAvg(emap.partFlux);
309     }
310    
311     portCnt++;
312 rschregle 2.12 } while (portCnt < numPhotonPorts);
313    
314 greg 2.1 if (srcFlux [srcIdx] < FTINY) {
315     sprintf(errmsg, "source %s has zero emission",
316     source [srcIdx].so -> oname);
317     error(WARNING, errmsg);
318     }
319 rschregle 2.12 }
320 rschregle 2.14
321     /* Allocate & init per-subprocess primary heap files */
322     primaryHeap = calloc(numProc, sizeof(FILE*));
323     primaryHeapFname = calloc(numProc, sizeof(char*));
324     primaryOfs = calloc(numProc, sizeof(PhotonPrimaryIdx));
325     if (!primaryHeap || !primaryHeapFname || !primaryOfs)
326     error(SYSTEM, "failed primary heap allocation in "
327     "distribPhotonContrib");
328    
329     for (proc = 0; proc < numProc; proc++) {
330     primaryHeapFname [proc] = malloc(PMAP_TMPFNLEN);
331     if (!primaryHeapFname [proc])
332     error(SYSTEM, "failed primary heap file allocation in "
333     "distribPhotonContrib");
334    
335     mktemp(strcpy(primaryHeapFname [proc], PMAP_TMPFNAME));
336     if (!(primaryHeap [proc] = fopen(primaryHeapFname [proc], "w+b")))
337     error(SYSTEM, "failed opening primary heap file in "
338     "distribPhotonContrib");
339     }
340 rschregle 2.12
341 rschregle 2.14 /* Record start time for progress reports */
342     repStartTime = time(NULL);
343 rschregle 2.12
344 rschregle 2.14 if (verbose) {
345     sprintf(errmsg, "\nPhoton distribution @ %d procs\n", numProc);
346     eputs(errmsg);
347     }
348 rschregle 2.12
349     /* MAIN LOOP */
350     for (proc = 0; proc < numProc; proc++) {
351 rschregle 2.14 #if NIX
352 rschregle 2.12 if (!(pid = fork())) {
353 rschregle 2.14 /* SUBPROCESS ENTERS HERE; opened and mmapped files inherited */
354     #else
355     if (1) {
356     /* No subprocess under Windoze */
357     #endif
358 rschregle 2.12 /* Local photon counters for this subprocess */
359     unsigned long lastNumPhotons = 0, localNumEmitted = 0;
360 rschregle 2.14 double photonFluxSum = 0; /* Accum. photon flux */
361 rschregle 2.12
362     /* Seed RNGs from PID for decorellated photon distribution */
363     pmapSeed(randSeed + proc, partState);
364 rschregle 2.16 pmapSeed(randSeed + (proc + 1) % numProc, emitState);
365     pmapSeed(randSeed + (proc + 2) % numProc, cntState);
366     pmapSeed(randSeed + (proc + 3) % numProc, mediumState);
367     pmapSeed(randSeed + (proc + 4) % numProc, scatterState);
368     pmapSeed(randSeed + (proc + 5) % numProc, rouletteState);
369    
370     #ifdef PMAP_SIGUSR
371     double partNumEmit;
372     unsigned long partEmitCnt;
373     double srcPhotonFlux, avgPhotonFlux;
374     unsigned portCnt, passCnt, prePassCnt;
375     float srcPreDistrib;
376     double srcNumEmit; /* # to emit from source */
377     unsigned long srcNumDistrib; /* # stored */
378    
379     void sigUsrDiags()
380     /* Loop diags via SIGUSR1 */
381     {
382     sprintf(errmsg,
383     "********************* Proc %d Diags *********************\n"
384     "srcIdx = %d (%s)\nportCnt = %d (%s)\npassCnt = %d\n"
385     "srcFlux = %f\nsrcPhotonFlux = %f\navgPhotonFlux = %f\n"
386     "partNumEmit = %f\npartEmitCnt = %lu\n\n",
387     proc, srcIdx, findmaterial(source [srcIdx].so) -> oname,
388     portCnt, photonPorts [portCnt].so -> oname,
389     passCnt, srcFlux [srcIdx], srcPhotonFlux, avgPhotonFlux,
390     partNumEmit, partEmitCnt);
391     eputs(errmsg);
392     fflush(stderr);
393     }
394     #endif
395    
396     #if PMAP_SIGUSR
397     signal(SIGUSR1, sigUsrDiags);
398     #endif
399 rschregle 2.12
400 rschregle 2.16 /* Output child process PID after random delay to prevent corrupted
401     * console output due to race condition */
402     usleep(1e6 * pmapRandom(rouletteState));
403     fprintf(stderr, "Proc %d: PID = %d\n", proc, getpid());
404     /* Allow time for debugger to attach to child process */
405     sleep(10);
406    
407 rschregle 2.12 /* =============================================================
408 greg 2.1 * 2-PASS PHOTON DISTRIBUTION
409     * Pass 1 (pre): emit fraction of target photon count
410 rschregle 2.12 * Pass 2 (main): based on outcome of pass 1, estimate remaining
411     * number of photons to emit to approximate target
412     * count
413     * ============================================================= */
414     for (srcIdx = 0; srcIdx < nsources; srcIdx++) {
415 rschregle 2.16 #ifndef PMAP_SIGUSR
416 rschregle 2.12 unsigned portCnt, passCnt = 0, prePassCnt = 0;
417     float srcPreDistrib = preDistrib;
418     double srcNumEmit = 0; /* # to emit from source */
419     unsigned long srcNumDistrib = pm -> numPhotons; /* # stored */
420 rschregle 2.16 #else
421     passCnt = prePassCnt = 0;
422     srcPreDistrib = preDistrib;
423     srcNumEmit = 0; /* # to emit from source */
424     srcNumDistrib = pm -> numPhotons; /* # stored */
425     #endif
426 rschregle 2.12
427     if (srcFlux [srcIdx] < FTINY)
428     continue;
429    
430     while (passCnt < 2) {
431     if (!passCnt) {
432     /* INIT PASS 1 */
433 rschregle 2.16 if (++prePassCnt > maxPreDistrib) {
434 rschregle 2.12 /* Warn if no photons contributed after sufficient
435 rschregle 2.14 * iterations; only output from subprocess 0 to reduce
436     * console clutter */
437 rschregle 2.16 if (!proc) {
438     sprintf(errmsg,
439     "source %s: too many prepasses, skipped",
440     source [srcIdx].so -> oname);
441     error(WARNING, errmsg);
442     }
443    
444 rschregle 2.12 break;
445     }
446    
447     /* Num to emit is fraction of target count */
448     srcNumEmit = srcPreDistrib * srcDistribTarget;
449 greg 2.1 }
450 rschregle 2.12 else {
451     /* INIT PASS 2 */
452 rschregle 2.16 #ifndef PMAP_SIGUSR
453 rschregle 2.12 double srcPhotonFlux, avgPhotonFlux;
454 rschregle 2.16 #endif
455 rschregle 2.12
456     /* Based on the outcome of the predistribution we can now
457     * figure out how many more photons we have to emit from
458     * the current source to meet the target count,
459     * srcDistribTarget. This value is clamped to 0 in case
460     * the target has already been exceeded in pass 1.
461     * srcNumEmit and srcNumDistrib is the number of photons
462     * emitted and distributed (stored) from the current
463     * source in pass 1, respectively. */
464     srcNumDistrib = pm -> numPhotons - srcNumDistrib;
465     srcNumEmit *= srcNumDistrib
466     ? max(srcDistribTarget/srcNumDistrib, 1) - 1
467     : 0;
468    
469     if (!srcNumEmit)
470     /* No photons left to distribute in main pass */
471     break;
472 greg 2.1
473 rschregle 2.12 srcPhotonFlux = srcFlux [srcIdx] / srcNumEmit;
474     avgPhotonFlux = photonFluxSum / (srcIdx + 1);
475    
476 rschregle 2.16 if (avgPhotonFlux > FTINY &&
477 rschregle 2.12 srcPhotonFlux / avgPhotonFlux < FTINY) {
478     /* Skip source if its photon flux is grossly below the
479 rschregle 2.14 * running average, indicating negligible contributions
480     * at the expense of excessive distribution time; only
481     * output from subproc 0 to reduce console clutter */
482 rschregle 2.16 if (!proc) {
483     sprintf(errmsg,
484     "source %s: itsy bitsy photon flux, skipped",
485     source [srcIdx].so -> oname);
486     error(WARNING, errmsg);
487     }
488    
489     srcNumEmit = 0; /* Or just break??? */
490 greg 2.1 }
491 rschregle 2.12
492     /* Update sum of photon flux per light source */
493     photonFluxSum += srcPhotonFlux;
494 greg 2.1 }
495 rschregle 2.14
496 rschregle 2.12 portCnt = 0;
497     do { /* Need at least one iteration if no ports! */
498     emap.src = source + srcIdx;
499     emap.port = emap.src -> sflags & SDISTANT
500     ? photonPorts + portCnt : NULL;
501     photonPartition [emap.src -> so -> otype] (&emap);
502 rschregle 2.14
503     if (verbose && !proc) {
504     /* Output from subproc 0 only to avoid race condition
505     * on console I/O */
506 rschregle 2.12 if (!passCnt)
507 rschregle 2.14 sprintf(errmsg, "\tPREPASS %d on source %s ",
508     prePassCnt, source [srcIdx].so -> oname);
509 rschregle 2.12 else
510 rschregle 2.14 sprintf(errmsg, "\tMAIN PASS on source %s ",
511     source [srcIdx].so -> oname);
512    
513 rschregle 2.12 if (emap.port) {
514     sprintf(errmsg2, "via port %s ",
515     photonPorts [portCnt].so -> oname);
516     strcat(errmsg, errmsg2);
517     }
518 rschregle 2.14
519 rschregle 2.12 sprintf(errmsg2, "(%lu partitions)\n",
520     emap.numPartitions);
521 rschregle 2.14 strcat(errmsg, errmsg2);
522 rschregle 2.12 eputs(errmsg);
523 rschregle 2.14 #if NIX
524 rschregle 2.12 fflush(stderr);
525 rschregle 2.14 #endif
526     }
527 greg 2.1
528 rschregle 2.12 for (emap.partitionCnt = 0; emap.partitionCnt < emap.numPartitions;
529     emap.partitionCnt++) {
530 rschregle 2.16 #ifndef PMAP_SIGUSR
531 rschregle 2.12 double partNumEmit;
532     unsigned long partEmitCnt;
533 rschregle 2.16 #endif
534 greg 2.1
535 rschregle 2.12 /* Get photon origin within current source partishunn
536     * and build emission map */
537     photonOrigin [emap.src -> so -> otype] (&emap);
538     initPhotonEmission(&emap, pdfSamples);
539    
540     /* Number of photons to emit from ziss partishunn;
541     * scale according to its normalised contribushunn to
542     * the emitted source flux */
543     partNumEmit = srcNumEmit * colorAvg(emap.partFlux) /
544     srcFlux [srcIdx];
545     partEmitCnt = (unsigned long)partNumEmit;
546    
547     /* Probabilistically account for fractional photons */
548     if (pmapRandom(cntState) < partNumEmit - partEmitCnt)
549     partEmitCnt++;
550    
551     /* Update local and shared global emission counter */
552 rschregle 2.14 photonCnt [PHOTONCNT_NUMEMIT(srcIdx)] += partEmitCnt;
553 rschregle 2.12 localNumEmitted += partEmitCnt;
554    
555 rschregle 2.14 /* Integer counter avoids FP rounding errors during
556     * iteration */
557 rschregle 2.12 while (partEmitCnt--) {
558     RAY photonRay;
559 greg 2.1
560 rschregle 2.12 /* Emit photon according to PDF (if any), allocate
561     * associated primary ray, and trace through scene
562     * until absorbed/leaked; emitPhoton() sets the
563     * emitting light source index in photonRay */
564     emitPhoton(&emap, &photonRay);
565 rschregle 2.14 #if 1
566     if (emap.port)
567     /* !!! PHOTON PORT REJECTION SAMPLING HACK: set
568     * !!! photon port as fake hit object for
569     * !!! primary ray to check for intersection in
570     * !!! tracePhoton() */
571     photonRay.ro = emap.port -> so;
572     #endif
573 rschregle 2.12 newPhotonPrimary(pm, &photonRay, primaryHeap[proc]);
574     /* Set subprocess index in photonRay for post-
575     * distrib primary index linearisation; this is
576     * propagated with the primary index in photonRay
577     * and set for photon hits by newPhoton() */
578     PMAP_SETRAYPROC(&photonRay, proc);
579     tracePhoton(&photonRay);
580     }
581 greg 2.1
582 rschregle 2.12 /* Update shared global photon count */
583     photonCnt [PHOTONCNT_NUMPHOT] += pm -> numPhotons -
584     lastNumPhotons;
585     lastNumPhotons = pm -> numPhotons;
586 rschregle 2.14 #if !NIX
587     /* Synchronous progress report on Windoze */
588     if (!proc && photonRepTime > 0 &&
589     time(NULL) >= repLastTime + photonRepTime) {
590     unsigned s;
591     repComplete = pm -> distribTarget * numProc;
592     repProgress = photonCnt [PHOTONCNT_NUMPHOT];
593    
594     for (repEmitted = 0, s = 0; s < nsources; s++)
595     repEmitted += photonCnt [PHOTONCNT_NUMEMIT(s)];
596    
597     pmapDistribReport();
598     }
599     #endif
600 greg 2.1 }
601 rschregle 2.12
602     portCnt++;
603     } while (portCnt < numPhotonPorts);
604    
605 rschregle 2.14 if (pm -> numPhotons == srcNumDistrib) {
606 rschregle 2.12 /* Double predistrib factor in case no photons were stored
607     * for this source and redo pass 1 */
608     srcPreDistrib *= 2;
609 rschregle 2.14 }
610 rschregle 2.12 else {
611     /* Now do pass 2 */
612     passCnt++;
613 greg 2.1 }
614     }
615 rschregle 2.12 }
616    
617     /* Flush heap buffa one final time to prevent data corruption */
618 rschregle 2.14 flushPhotonHeap(pm);
619 rschregle 2.12 /* Flush final photon primary to primary heap file */
620     newPhotonPrimary(pm, NULL, primaryHeap [proc]);
621 rschregle 2.14 /* Heap files closed automatically on exit
622     fclose(pm -> heap);
623     fclose(primaryHeap [proc]); */
624 rschregle 2.12
625     #ifdef DEBUG_PMAP
626 rschregle 2.14 sprintf(errmsg, "Proc %d total %ld photons\n", proc,
627 rschregle 2.12 pm -> numPhotons);
628     eputs(errmsg);
629 rschregle 2.14 fflush(stderr);
630 rschregle 2.12 #endif
631    
632 rschregle 2.16 #ifdef PMAP_SIGUSR
633     signal(SIGUSR1, SIG_DFL);
634     #endif
635    
636 rschregle 2.14 #if NIX
637     /* Terminate subprocess */
638 rschregle 2.12 exit(0);
639 rschregle 2.14 #endif
640 greg 2.1 }
641 rschregle 2.12 else if (pid < 0)
642     error(SYSTEM, "failed to fork subprocess in distribPhotonContrib");
643     }
644    
645 rschregle 2.14 #if NIX
646 rschregle 2.12 /* PARENT PROCESS CONTINUES HERE */
647     #ifdef SIGCONT
648 rschregle 2.14 /* Enable progress report signal handler */
649 rschregle 2.12 signal(SIGCONT, pmapDistribReport);
650     #endif
651     /* Wait for subprocesses to complete while reporting progress */
652     proc = numProc;
653     while (proc) {
654     while (waitpid(-1, &stat, WNOHANG) > 0) {
655     /* Subprocess exited; check status */
656     if (!WIFEXITED(stat) || WEXITSTATUS(stat))
657     error(USER, "failed photon distribution");
658    
659     --proc;
660     }
661    
662     /* Nod off for a bit and update progress */
663     sleep(1);
664 rschregle 2.14
665     /* Asynchronous progress report from shared subprocess counters */
666     repComplete = pm -> distribTarget * numProc;
667     repProgress = photonCnt [PHOTONCNT_NUMPHOT];
668 rschregle 2.12
669     for (repEmitted = 0, srcIdx = 0; srcIdx < nsources; srcIdx++)
670     repEmitted += photonCnt [PHOTONCNT_NUMEMIT(srcIdx)];
671    
672     /* Get global photon count from shmem updated by subprocs */
673     pm -> numPhotons = photonCnt [PHOTONCNT_NUMPHOT];
674    
675     if (photonRepTime > 0 && time(NULL) >= repLastTime + photonRepTime)
676     pmapDistribReport();
677     #ifdef SIGCONT
678     else signal(SIGCONT, pmapDistribReport);
679     #endif
680 greg 2.1 }
681 rschregle 2.14 #endif /* NIX */
682 greg 2.1
683     /* ================================================================
684     * POST-DISTRIBUTION - Set photon flux and build kd-tree, etc.
685     * ================================================================ */
686 rschregle 2.12 #ifdef SIGCONT
687 rschregle 2.14 /* Reset signal handler */
688 rschregle 2.12 signal(SIGCONT, SIG_DFL);
689     #endif
690 greg 2.1 free(emap.samples);
691    
692 rschregle 2.12 if (!pm -> numPhotons)
693 rschregle 2.14 error(USER, "empty contribution photon map");
694 greg 2.1
695 rschregle 2.12 /* Load per-subprocess primary rays into pm -> primary array */
696 rschregle 2.14 /* Dumb compilers apparently need the char** cast */
697     pm -> numPrimary = buildPrimaries(pm, primaryHeap,
698     (char**)primaryHeapFname,
699     primaryOfs, numProc);
700 rschregle 2.12 if (!pm -> numPrimary)
701 greg 2.1 error(INTERNAL, "no primary rays in contribution photon map");
702 rschregle 2.12
703     /* Set photon flux per source */
704     for (srcIdx = 0; srcIdx < nsources; srcIdx++)
705     srcFlux [srcIdx] /= photonCnt [PHOTONCNT_NUMEMIT(srcIdx)];
706 rschregle 2.14 #if NIX
707 rschregle 2.12 /* Photon counters no longer needed, unmap shared memory */
708     munmap(photonCnt, sizeof(*photonCnt));
709     close(shmFile);
710 rschregle 2.14 unlink(shmFname);
711 rschregle 2.12 #else
712 rschregle 2.14 free(photonCnt);
713 rschregle 2.12 #endif
714 greg 2.1
715 rschregle 2.14 if (verbose) {
716     eputs("\nBuilding contribution photon map...\n");
717     #if NIX
718 greg 2.1 fflush(stderr);
719 rschregle 2.14 #endif
720 greg 2.1 }
721 rschregle 2.12
722     /* Build underlying data structure; heap is destroyed */
723 rschregle 2.14 buildPhotonMap(pm, srcFlux, primaryOfs, numProc);
724    
725     /* Free per-subprocess primary heap files */
726     for (proc = 0; proc < numProc; proc++)
727     free(primaryHeapFname [proc]);
728    
729     free(primaryHeapFname);
730     free(primaryHeap);
731     free(primaryOfs);
732    
733     if (verbose)
734     eputs("\n");
735 rschregle 2.12 }