ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/rt/source.c
Revision: 2.39
Committed: Wed Dec 31 02:03:08 2003 UTC (20 years, 4 months ago) by greg
Content type: text/plain
Branch: MAIN
Changes since 2.38: +3 -1 lines
Log Message:
Minor fixes to source occlusion cache

File Contents

# Content
1 #ifndef lint
2 static const char RCSid[] = "$Id: source.c,v 2.38 2003/12/31 01:50:02 greg Exp $";
3 #endif
4 /*
5 * source.c - routines dealing with illumination sources.
6 *
7 * External symbols declared in source.h
8 */
9
10 #include "copyright.h"
11
12 #include "ray.h"
13
14 #include "otypes.h"
15
16 #include "otspecial.h"
17
18 #include "source.h"
19
20 #include "random.h"
21
22 extern double ssampdist; /* scatter sampling distance */
23
24 #ifndef MAXSSAMP
25 #define MAXSSAMP 16 /* maximum samples per ray */
26 #endif
27
28 /*
29 * Structures used by direct()
30 */
31
32 typedef struct {
33 int sno; /* source number */
34 FVECT dir; /* source direction */
35 COLOR coef; /* material coefficient */
36 COLOR val; /* contribution */
37 } CONTRIB; /* direct contribution */
38
39 typedef struct {
40 int sndx; /* source index (to CONTRIB array) */
41 float brt; /* brightness (for comparison) */
42 } CNTPTR; /* contribution pointer */
43
44 static CONTRIB *srccnt; /* source contributions in direct() */
45 static CNTPTR *cntord; /* source ordering in direct() */
46 static int maxcntr = 0; /* size of contribution arrays */
47
48
49 OBJREC * /* find an object's actual material */
50 findmaterial(register OBJREC *o)
51 {
52 while (!ismaterial(o->otype)) {
53 if (ismixture(o->otype))
54 return(NULL); /* reject mixed materials */
55 if (o->otype == MOD_ALIAS && o->oargs.nsargs) {
56 OBJECT aobj;
57 OBJREC *ao;
58 aobj = lastmod(objndx(o), o->oargs.sarg[0]);
59 if (aobj < 0)
60 objerror(o, USER, "bad reference");
61 ao = objptr(aobj);
62 if (ismaterial(ao->otype))
63 return(ao);
64 }
65 if (o->omod == OVOID)
66 return(NULL);
67 o = objptr(o->omod);
68 }
69 return(o);
70 }
71
72
73 void
74 marksources() /* find and mark source objects */
75 {
76 int foundsource = 0;
77 int i;
78 register OBJREC *o, *m;
79 register int ns;
80 /* initialize dispatch table */
81 initstypes();
82 /* find direct sources */
83 for (i = 0; i < nsceneobjs; i++) {
84
85 o = objptr(i);
86
87 if (!issurface(o->otype) || o->omod == OVOID)
88 continue;
89 /* find material */
90 m = findmaterial(o);
91 if (m == NULL || !islight(m->otype))
92 continue; /* not source modifier */
93
94 if (m->oargs.nfargs != (m->otype == MAT_GLOW ? 4 :
95 m->otype == MAT_SPOT ? 7 : 3))
96 objerror(m, USER, "bad # arguments");
97
98 if (m->otype == MAT_GLOW &&
99 o->otype != OBJ_SOURCE &&
100 m->oargs.farg[3] <= FTINY)
101 continue; /* don't bother */
102 if (m->oargs.farg[0] <= FTINY && m->oargs.farg[1] <= FTINY &&
103 m->oargs.farg[2] <= FTINY)
104 continue; /* don't bother */
105
106 if (sfun[o->otype].of == NULL ||
107 sfun[o->otype].of->setsrc == NULL)
108 objerror(o, USER, "illegal material");
109
110 if ((ns = newsource()) < 0)
111 goto memerr;
112
113 setsource(&source[ns], o);
114
115 if (m->otype == MAT_GLOW) {
116 source[ns].sflags |= SPROX;
117 source[ns].sl.prox = m->oargs.farg[3];
118 if (source[ns].sflags & SDISTANT)
119 source[ns].sflags |= SSKIP;
120 } else if (m->otype == MAT_SPOT) {
121 source[ns].sflags |= SSPOT;
122 if ((source[ns].sl.s = makespot(m)) == NULL)
123 goto memerr;
124 if (source[ns].sflags & SFLAT &&
125 !checkspot(source[ns].sl.s,source[ns].snorm)) {
126 objerror(o, WARNING,
127 "invalid spotlight direction");
128 source[ns].sflags |= SSKIP;
129 }
130 }
131 #if SHADCACHE
132 source[ns].obscache = NULL;
133 #endif
134 if (!(source[ns].sflags & SSKIP))
135 foundsource++;
136 }
137 if (!foundsource) {
138 error(WARNING, "no light sources found");
139 return;
140 }
141 markvirtuals(); /* find and add virtual sources */
142 /* allocate our contribution arrays */
143 maxcntr = nsources + MAXSPART; /* start with this many */
144 srccnt = (CONTRIB *)malloc(maxcntr*sizeof(CONTRIB));
145 cntord = (CNTPTR *)malloc(maxcntr*sizeof(CNTPTR));
146 if ((srccnt == NULL) | (cntord == NULL))
147 goto memerr;
148 return;
149 memerr:
150 error(SYSTEM, "out of memory in marksources");
151 }
152
153
154 void
155 freesources() /* free all source structures */
156 {
157 if (nsources > 0) {
158 #if SHADCACHE
159 while (nsources--)
160 freeobscache(&source[nsources]);
161 #endif
162 free((void *)source);
163 source = NULL;
164 nsources = 0;
165 }
166 if (maxcntr <= 0)
167 return;
168 free((void *)srccnt);
169 srccnt = NULL;
170 free((void *)cntord);
171 cntord = NULL;
172 maxcntr = 0;
173 }
174
175
176 int
177 srcray( /* send a ray to a source, return domega */
178 register RAY *sr, /* returned source ray */
179 RAY *r, /* ray which hit object */
180 SRCINDEX *si /* source sample index */
181 )
182 {
183 double d; /* distance to source */
184 register SRCREC *srcp;
185
186 rayorigin(sr, r, SHADOW, 1.0); /* ignore limits */
187
188 while ((d = nextssamp(sr, si)) != 0.0) {
189 sr->rsrc = si->sn; /* remember source */
190 srcp = source + si->sn;
191 if (srcp->sflags & SDISTANT) {
192 if (srcp->sflags & SSPOT && spotout(sr, srcp->sl.s))
193 continue;
194 return(1); /* sample OK */
195 }
196 /* local source */
197 /* check proximity */
198 if (srcp->sflags & SPROX && d > srcp->sl.prox)
199 continue;
200 /* check angle */
201 if (srcp->sflags & SSPOT) {
202 if (spotout(sr, srcp->sl.s))
203 continue;
204 /* adjust solid angle */
205 si->dom *= d*d;
206 d += srcp->sl.s->flen;
207 si->dom /= d*d;
208 }
209 return(1); /* sample OK */
210 }
211 return(0); /* no more samples */
212 }
213
214
215 void
216 srcvalue( /* punch ray to source and compute value */
217 register RAY *r
218 )
219 {
220 register SRCREC *sp;
221
222 sp = &source[r->rsrc];
223 if (sp->sflags & SVIRTUAL) { /* virtual source */
224 /* check intersection */
225 if (!(*ofun[sp->so->otype].funp)(sp->so, r))
226 return;
227 if (!rayshade(r, r->ro->omod)) /* compute contribution */
228 goto nomat;
229 rayparticipate(r);
230 return;
231 }
232 /* compute intersection */
233 if (sp->sflags & SDISTANT ? sourcehit(r) :
234 (*ofun[sp->so->otype].funp)(sp->so, r)) {
235 if (sp->sa.success >= 0)
236 sp->sa.success++;
237 if (!rayshade(r, r->ro->omod)) /* compute contribution */
238 goto nomat;
239 rayparticipate(r);
240 return;
241 }
242 /* we missed our mark! */
243 if (sp->sa.success < 0)
244 return; /* bitched already */
245 sp->sa.success -= AIMREQT;
246 if (sp->sa.success >= 0)
247 return; /* leniency */
248 sprintf(errmsg, "aiming failure for light source \"%s\"",
249 sp->so->oname);
250 error(WARNING, errmsg); /* issue warning */
251 return;
252 nomat:
253 objerror(r->ro, USER, "material not found");
254 }
255
256
257 int
258 sourcehit( /* check to see if ray hit distant source */
259 register RAY *r
260 )
261 {
262 int first, last;
263 register int i;
264
265 if (r->rsrc >= 0) { /* check only one if aimed */
266 first = last = r->rsrc;
267 } else { /* otherwise check all */
268 first = 0; last = nsources-1;
269 }
270 for (i = first; i <= last; i++)
271 if ((source[i].sflags & (SDISTANT|SVIRTUAL)) == SDISTANT)
272 /*
273 * Check to see if ray is within
274 * solid angle of source.
275 */
276 if (2.0*PI * (1.0 - DOT(source[i].sloc,r->rdir))
277 <= source[i].ss2) {
278 r->ro = source[i].so;
279 if (!(source[i].sflags & SSKIP))
280 break;
281 }
282
283 if (r->ro != NULL) {
284 r->robj = objndx(r->ro);
285 for (i = 0; i < 3; i++)
286 r->ron[i] = -r->rdir[i];
287 r->rod = 1.0;
288 r->pert[0] = r->pert[1] = r->pert[2] = 0.0;
289 r->uv[0] = r->uv[1] = 0.0;
290 r->rox = NULL;
291 return(1);
292 }
293 return(0);
294 }
295
296
297 #if SHADCACHE /* preemptive shadow checking */
298 #define ABS(x) ((x)>0 ? (x) : -(x))
299
300 static void /* find closest blockers to source */
301 initobscache(SRCREC *srcp)
302 {
303 int i;
304 int cachelen;
305
306 if (srcp->sflags & SDISTANT)
307 cachelen = 4*SHADCACHE*SHADCACHE;
308 else if (srcp->sflags & SFLAT)
309 cachelen = SHADCACHE*SHADCACHE*3 + (SHADCACHE&1)*SHADCACHE*4;
310 else /* spherical distribution */
311 cachelen = SHADCACHE*SHADCACHE*6;
312 /* allocate cache */
313 DCHECK(srcp->obscache != NULL,
314 CONSISTENCY, "initobscache() called twice");
315 srcp->obscache = (OBSCACHE *)malloc(sizeof(OBSCACHE) +
316 sizeof(OBJECT)*(cachelen-1));
317 if (srcp->obscache == NULL)
318 error(SYSTEM, "out of memory in initobscache()");
319 /* set parameters */
320 if (srcp->sflags & SDISTANT) {
321 int ax, ax1, ax2;
322 RREAL amax = 0;
323 for (ax1 = 3; ax1--; )
324 if (ABS(srcp->sloc[ax1]) > amax) {
325 amax = ABS(srcp->sloc[ax1]);
326 ax = ax1;
327 }
328 srcp->obscache->p.d.ax = ax;
329 ax1 = (ax+1)%3;
330 ax2 = (ax+2)%3;
331 VCOPY(srcp->obscache->p.d.o, thescene.cuorg);
332 if (srcp->sloc[ax] > 0)
333 srcp->obscache->p.d.o[ax] += thescene.cusize;
334 if (srcp->sloc[ax1] < 0)
335 srcp->obscache->p.d.o[ax1] += thescene.cusize *
336 srcp->sloc[ax1] / ABS(srcp->sloc[ax]);
337 if (srcp->sloc[ax2] < 0)
338 srcp->obscache->p.d.o[ax2] += thescene.cusize *
339 srcp->sloc[ax2] / ABS(srcp->sloc[ax]);
340 srcp->obscache->p.d.e1 = (1.-FTINY) / (thescene.cusize*(1. +
341 fabs(srcp->sloc[ax1]/srcp->sloc[ax])));
342 srcp->obscache->p.d.e2 = (1.-FTINY) / (thescene.cusize*(1. +
343 fabs(srcp->sloc[ax2]/srcp->sloc[ax])));
344 } else if (srcp->sflags & SFLAT) {
345 VCOPY(srcp->obscache->p.f.u, srcp->ss[SU]);
346 normalize(srcp->obscache->p.f.u);
347 fcross(srcp->obscache->p.f.v,
348 srcp->snorm, srcp->obscache->p.f.u);
349 }
350 /* XXX Should cast rays from source */
351 for (i = cachelen; i--; )
352 srcp->obscache->obs[i] = OVOID;
353 }
354
355
356 static OBJECT * /* return occluder cache entry */
357 srcobstructp(register RAY *r)
358 {
359 static OBJECT noobs;
360 SRCREC *srcp;
361 int ondx;
362
363 DCHECK(r->rsrc < 0, CONSISTENCY,
364 "srcobstructp() called with unaimed ray");
365 noobs = OVOID;
366 srcp = &source[r->rsrc];
367 if (srcp->obscache == NULL) /* initialize cache */
368 initobscache(srcp);
369 /* compute cache index */
370 if (srcp->sflags & SDISTANT) {
371 int ax, ax1, ax2;
372 double t;
373 ax = srcp->obscache->p.d.ax;
374 if ((ax1 = ax+1) >= 3) ax1 -= 3;
375 if ((ax2 = ax+2) >= 3) ax2 -= 3;
376 t = (srcp->obscache->p.d.o[ax] - r->rorg[ax]) / srcp->sloc[ax];
377 if (t <= FTINY)
378 return &noobs; /* could happen if ray is outside */
379 ondx = 2*SHADCACHE*(int)(2*SHADCACHE*srcp->obscache->p.d.e1 *
380 (r->rorg[ax1] + t*srcp->sloc[ax1] -
381 srcp->obscache->p.d.o[ax1]));
382 ondx += (int)(2*SHADCACHE*srcp->obscache->p.d.e2 *
383 (r->rorg[ax2] + t*srcp->sloc[ax2] -
384 srcp->obscache->p.d.o[ax2]));
385 if (ondx < 0 | ondx >= 4*SHADCACHE*SHADCACHE)
386 return &nobs; /* could happen if ray is outside */
387 } else if (srcp->sflags & SFLAT) {
388 FVECT sd;
389 RREAL sd0m, sd1m;
390 sd[0] = -DOT(r->rdir, srcp->obscache->p.f.u);
391 sd[1] = -DOT(r->rdir, srcp->obscache->p.f.v);
392 sd[2] = -DOT(r->rdir, srcp->snorm);
393 if (sd[2] < 0)
394 return &noobs; /* shouldn't happen */
395 sd0m = ABS(sd[0]);
396 sd1m = ABS(sd[1]);
397 if (sd[2] >= sd0m && sd[2] >= sd1m) {
398 ondx = SHADCACHE*(int)(SHADCACHE*(.5-FTINY) *
399 (1. + sd[0]/sd[2]));
400 ondx += (int)(SHADCACHE*(.5-FTINY) *
401 (1. + sd[1]/sd[2]));
402 } else if (sd0m >= sd1m) {
403 ondx = SHADCACHE*SHADCACHE;
404 if (sd[0] < 0)
405 ondx += ((SHADCACHE+1)>>1)*SHADCACHE;
406 ondx += SHADCACHE*(int)(SHADCACHE*(.5-FTINY) *
407 sd[2]/sd0m);
408 ondx += (int)(SHADCACHE*(.5-FTINY) *
409 (1. + sd[1]/sd0m));
410 } else /* sd1m > sd0m */ {
411 ondx = SHADCACHE*SHADCACHE +
412 ((SHADCACHE+1)>>1)*SHADCACHE*2;
413 if (sd[1] < 0)
414 ondx += ((SHADCACHE+1)>>1)*SHADCACHE;
415 ondx += SHADCACHE*(int)(SHADCACHE*(.5-FTINY) *
416 sd[2]/sd1m);
417 ondx += (int)(SHADCACHE*(.5-FTINY) *
418 (1. + sd[0]/sd1m));
419 }
420 } else /* spherical distribution */ {
421 int ax, ax1, ax2;
422 RREAL amax = 0;
423 for (ax1 = 3; ax1--; )
424 if (ABS(r->rdir[ax1]) > amax) {
425 amax = ABS(r->rdir[ax1]);
426 ax = ax1;
427 }
428 if ((ax1 = ax+1) >= 3) ax1 -= 3;
429 if ((ax2 = ax+2) >= 3) ax2 -= 3;
430 ondx = 2*SHADCACHE*SHADCACHE * ax;
431 if (r->rdir[ax] < 0)
432 ondx += SHADCACHE*SHADCACHE;
433 ondx += SHADCACHE*(int)(SHADCACHE*(.5-FTINY) *
434 (1. + r->rdir[ax1]/amax));
435 ondx += (int)(SHADCACHE*(.5-FTINY) *
436 (1. + r->rdir[ax2]/amax));
437 }
438 /* return cache pointer */
439 return(&srcp->obscache->obs[ondx]);
440 }
441
442
443 void /* free obstruction cache */
444 freeobscache(SRCREC *srcp)
445 {
446 if (srcp->obscache == NULL)
447 return;
448 free((void *)srcp->obscache);
449 srcp->obscache = NULL;
450 }
451
452
453 void /* record a source blocker */
454 srcblocker(register RAY *r)
455 {
456 OBJREC *m;
457
458 if (r->robj == OVOID || objptr(r->robj) != r->ro ||
459 isvolume(r->ro->otype))
460 return; /* don't record complex blockers */
461 m = findmaterial(r->ro);
462 if (m == NULL)
463 return; /* no material?! */
464 if (!(ofun[m->otype].flags & T_OPAQUE))
465 return; /* material not a reliable blocker */
466
467 *srcobstructp(r) = r->robj; /* else record obstructor */
468 }
469
470
471 int /* check ray against cached blocker */
472 srcblocked(RAY *r)
473 {
474 OBJECT obs = *srcobstructp(r);
475 OBJREC *op;
476
477 if (obs == OVOID)
478 return(0);
479 op = objptr(obs); /* check for intersection */
480 return ((*ofun[op->otype].funp)(op, r));
481 }
482
483 #endif
484
485
486 static int
487 cntcmp( /* contribution compare (descending) */
488 const void *p1,
489 const void *p2
490 )
491 {
492 register const CNTPTR *sc1 = (const CNTPTR *)p1;
493 register const CNTPTR *sc2 = (const CNTPTR *)p2;
494
495 if (sc1->brt > sc2->brt)
496 return(-1);
497 if (sc1->brt < sc2->brt)
498 return(1);
499 return(0);
500 }
501
502
503 void
504 direct( /* add direct component */
505 RAY *r, /* ray that hit surface */
506 void (*f)(), /* direct component coefficient function */
507 char *p /* data for f */
508 )
509 {
510 extern void (*trace)();
511 register int sn;
512 register CONTRIB *scp;
513 SRCINDEX si;
514 int nshadcheck, ncnts;
515 int nhits;
516 double prob, ourthresh, hwt;
517 RAY sr;
518 /* NOTE: srccnt and cntord global so no recursion */
519 if (nsources <= 0)
520 return; /* no sources?! */
521 /* potential contributions */
522 initsrcindex(&si);
523 for (sn = 0; srcray(&sr, r, &si); sn++) {
524 if (sn >= maxcntr) {
525 maxcntr = sn + MAXSPART;
526 srccnt = (CONTRIB *)realloc((void *)srccnt,
527 maxcntr*sizeof(CONTRIB));
528 cntord = (CNTPTR *)realloc((void *)cntord,
529 maxcntr*sizeof(CNTPTR));
530 if ((srccnt == NULL) | (cntord == NULL))
531 error(SYSTEM, "out of memory in direct");
532 }
533 cntord[sn].sndx = sn;
534 scp = srccnt + sn;
535 scp->sno = sr.rsrc;
536 /* compute coefficient */
537 (*f)(scp->coef, p, sr.rdir, si.dom);
538 cntord[sn].brt = bright(scp->coef);
539 if (cntord[sn].brt <= 0.0)
540 continue;
541 #if SHADCACHE
542 /* check shadow cache */
543 if (si.np == 1 && srcblocked(&sr)) {
544 cntord[sn].brt = 0.0;
545 continue;
546 }
547 #endif
548 VCOPY(scp->dir, sr.rdir);
549 /* compute potential */
550 sr.revf = srcvalue;
551 rayvalue(&sr);
552 copycolor(scp->val, sr.rcol);
553 multcolor(scp->val, scp->coef);
554 cntord[sn].brt = bright(scp->val);
555 }
556 /* sort contributions */
557 qsort(cntord, sn, sizeof(CNTPTR), cntcmp);
558 { /* find last */
559 register int l, m;
560
561 ncnts = l = sn;
562 sn = 0;
563 while ((m = (sn + ncnts) >> 1) != l) {
564 if (cntord[m].brt > 0.0)
565 sn = m;
566 else
567 ncnts = m;
568 l = m;
569 }
570 }
571 if (ncnts == 0)
572 return; /* no contributions! */
573 /* accumulate tail */
574 for (sn = ncnts-1; sn > 0; sn--)
575 cntord[sn-1].brt += cntord[sn].brt;
576 /* compute number to check */
577 nshadcheck = pow((double)ncnts, shadcert) + .5;
578 /* modify threshold */
579 ourthresh = shadthresh / r->rweight;
580 /* test for shadows */
581 for (nhits = 0, hwt = 0.0, sn = 0; sn < ncnts;
582 hwt += (double)source[scp->sno].nhits /
583 (double)source[scp->sno].ntests,
584 sn++) {
585 /* check threshold */
586 if ((sn+nshadcheck>=ncnts ? cntord[sn].brt :
587 cntord[sn].brt-cntord[sn+nshadcheck].brt)
588 < ourthresh*bright(r->rcol))
589 break;
590 scp = srccnt + cntord[sn].sndx;
591 /* test for hit */
592 rayorigin(&sr, r, SHADOW, 1.0);
593 VCOPY(sr.rdir, scp->dir);
594 sr.rsrc = scp->sno;
595 /* keep statistics */
596 if (source[scp->sno].ntests++ > 0xfffffff0) {
597 source[scp->sno].ntests >>= 1;
598 source[scp->sno].nhits >>= 1;
599 }
600 if (localhit(&sr, &thescene) &&
601 ( sr.ro != source[scp->sno].so ||
602 source[scp->sno].sflags & SFOLLOW )) {
603 /* follow entire path */
604 raycont(&sr);
605 rayparticipate(&sr);
606 if (trace != NULL)
607 (*trace)(&sr); /* trace execution */
608 if (bright(sr.rcol) <= FTINY) {
609 #if SHADCACHE
610 if ((scp <= srccnt || scp[-1].sno != scp->sno)
611 && (scp >= srccnt+ncnts ||
612 scp[1].sno != scp->sno))
613 srcblocker(&sr);
614 #endif
615 continue; /* missed! */
616 }
617 copycolor(scp->val, sr.rcol);
618 multcolor(scp->val, scp->coef);
619 }
620 /* add contribution if hit */
621 addcolor(r->rcol, scp->val);
622 nhits++;
623 source[scp->sno].nhits++;
624 }
625 /* source hit rate */
626 if (hwt > FTINY)
627 hwt = (double)nhits / hwt;
628 else
629 hwt = 0.5;
630 #ifdef DEBUG
631 sprintf(errmsg, "%d tested, %d untested, %f conditional hit rate\n",
632 sn, ncnts-sn, hwt);
633 eputs(errmsg);
634 #endif
635 /* add in untested sources */
636 for ( ; sn < ncnts; sn++) {
637 scp = srccnt + cntord[sn].sndx;
638 prob = hwt * (double)source[scp->sno].nhits /
639 (double)source[scp->sno].ntests;
640 if (prob > 1.0)
641 prob = 1.0;
642 scalecolor(scp->val, prob);
643 addcolor(r->rcol, scp->val);
644 }
645 }
646
647
648 void
649 srcscatter( /* compute source scattering into ray */
650 register RAY *r
651 )
652 {
653 int oldsampndx;
654 int nsamps;
655 RAY sr;
656 SRCINDEX si;
657 double t, d;
658 double re, ge, be;
659 COLOR cvext;
660 int i, j;
661
662 if (r->slights == NULL || r->slights[0] == 0
663 || r->gecc >= 1.-FTINY || r->rot >= FHUGE)
664 return;
665 if (ssampdist <= FTINY || (nsamps = r->rot/ssampdist + .5) < 1)
666 nsamps = 1;
667 #if MAXSSAMP
668 else if (nsamps > MAXSSAMP)
669 nsamps = MAXSSAMP;
670 #endif
671 oldsampndx = samplendx;
672 samplendx = random()&0x7fff; /* randomize */
673 for (i = r->slights[0]; i > 0; i--) { /* for each source */
674 for (j = 0; j < nsamps; j++) { /* for each sample position */
675 samplendx++;
676 t = r->rot * (j+frandom())/nsamps;
677 /* extinction */
678 re = t*colval(r->cext,RED);
679 ge = t*colval(r->cext,GRN);
680 be = t*colval(r->cext,BLU);
681 setcolor(cvext, re > 92. ? 0. : exp(-re),
682 ge > 92. ? 0. : exp(-ge),
683 be > 92. ? 0. : exp(-be));
684 if (intens(cvext) <= FTINY)
685 break; /* too far away */
686 sr.rorg[0] = r->rorg[0] + r->rdir[0]*t;
687 sr.rorg[1] = r->rorg[1] + r->rdir[1]*t;
688 sr.rorg[2] = r->rorg[2] + r->rdir[2]*t;
689 sr.rmax = 0.;
690 initsrcindex(&si); /* sample ray to this source */
691 si.sn = r->slights[i];
692 nopart(&si, &sr);
693 if (!srcray(&sr, NULL, &si) ||
694 sr.rsrc != r->slights[i])
695 continue; /* no path */
696 copycolor(sr.cext, r->cext);
697 copycolor(sr.albedo, r->albedo);
698 sr.gecc = r->gecc;
699 sr.slights = r->slights;
700 rayvalue(&sr); /* eval. source ray */
701 if (bright(sr.rcol) <= FTINY)
702 continue;
703 if (r->gecc <= FTINY) /* compute P(theta) */
704 d = 1.;
705 else {
706 d = DOT(r->rdir, sr.rdir);
707 d = 1. + r->gecc*r->gecc - 2.*r->gecc*d;
708 d = (1. - r->gecc*r->gecc) / (d*sqrt(d));
709 }
710 /* other factors */
711 d *= si.dom * r->rot / (4.*PI*nsamps);
712 multcolor(sr.rcol, r->cext);
713 multcolor(sr.rcol, r->albedo);
714 scalecolor(sr.rcol, d);
715 multcolor(sr.rcol, cvext);
716 addcolor(r->rcol, sr.rcol); /* add it in */
717 }
718 }
719 samplendx = oldsampndx;
720 }
721
722
723 /****************************************************************
724 * The following macros were separated from the m_light() routine
725 * because they are very nasty and difficult to understand.
726 */
727
728 /* illumblock *
729 *
730 * We cannot allow an illum to pass to another illum, because that
731 * would almost certainly constitute overcounting.
732 * However, we do allow an illum to pass to another illum
733 * that is actually going to relay to a virtual light source.
734 * We also prevent an illum from passing to a glow; this provides a
735 * convenient mechanism for defining detailed light source
736 * geometry behind (or inside) an effective radiator.
737 */
738
739 static int
740 weaksrcmat(int obj) /* identify material */
741 {
742 register OBJREC *o = objptr(obj);
743
744 while (!ismaterial(o->otype)) /* find material */
745 o = objptr(o->omod);
746 return((o->otype==MAT_ILLUM)|(o->otype==MAT_GLOW));
747 }
748
749 #define illumblock(m, r) (!(source[r->rsrc].sflags&SVIRTUAL) && \
750 r->rod > 0.0 && \
751 weaksrcmat(source[r->rsrc].so->omod))
752
753 /* wrongsource *
754 *
755 * This source is the wrong source (ie. overcounted) if we are
756 * aimed to a different source than the one we hit and the one
757 * we hit is not an illum that should be passed.
758 */
759
760 #define wrongsource(m, r) (r->rsrc>=0 && source[r->rsrc].so!=r->ro && \
761 (m->otype!=MAT_ILLUM || illumblock(m,r)))
762
763 /* distglow *
764 *
765 * A distant glow is an object that sometimes acts as a light source,
766 * but is too far away from the test point to be one in this case.
767 * (Glows with negative radii should NEVER participate in illumination.)
768 */
769
770 #define distglow(m, r, d) (m->otype==MAT_GLOW && \
771 m->oargs.farg[3] >= -FTINY && \
772 d > m->oargs.farg[3])
773
774 /* badcomponent *
775 *
776 * We must avoid counting light sources in the ambient calculation,
777 * since the direct component is handled separately. Therefore, any
778 * ambient ray which hits an active light source must be discarded.
779 * The same is true for stray specular samples, since the specular
780 * contribution from light sources is calculated separately.
781 */
782
783 #define badcomponent(m, r) (r->crtype&(AMBIENT|SPECULAR) && \
784 !(r->crtype&SHADOW || r->rod < 0.0 || \
785 /* not 100% correct */ distglow(m, r, r->rot)))
786
787 /* passillum *
788 *
789 * An illum passes to another material type when we didn't hit it
790 * on purpose (as part of a direct calculation), or it is relaying
791 * a virtual light source.
792 */
793
794 #define passillum(m, r) (m->otype==MAT_ILLUM && \
795 (r->rsrc<0 || source[r->rsrc].so!=r->ro || \
796 source[r->rsrc].sflags&SVIRTUAL))
797
798 /* srcignore *
799 *
800 * The -dv flag is normally on for sources to be visible.
801 */
802
803 #define srcignore(m, r) !(directvis || r->crtype&SHADOW || \
804 distglow(m, r, raydist(r,PRIMARY)))
805
806
807 int
808 m_light( /* ray hit a light source */
809 register OBJREC *m,
810 register RAY *r
811 )
812 {
813 /* check for over-counting */
814 if (badcomponent(m, r))
815 return(1);
816 if (wrongsource(m, r))
817 return(1);
818 /* check for passed illum */
819 if (passillum(m, r)) {
820 if (m->oargs.nsargs && strcmp(m->oargs.sarg[0], VOIDID))
821 return(rayshade(r,lastmod(objndx(m),m->oargs.sarg[0])));
822 raytrans(r);
823 return(1);
824 }
825 /* otherwise treat as source */
826 /* check for behind */
827 if (r->rod < 0.0)
828 return(1);
829 /* check for invisibility */
830 if (srcignore(m, r))
831 return(1);
832 /* check for outside spot */
833 if (m->otype==MAT_SPOT && spotout(r, makespot(m)))
834 return(1);
835 /* get distribution pattern */
836 raytexture(r, m->omod);
837 /* get source color */
838 setcolor(r->rcol, m->oargs.farg[0],
839 m->oargs.farg[1],
840 m->oargs.farg[2]);
841 /* modify value */
842 multcolor(r->rcol, r->pcol);
843 return(1);
844 }