ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/rt/source.c
Revision: 2.50
Committed: Tue Apr 19 01:15:07 2005 UTC (19 years ago) by greg
Content type: text/plain
Branch: MAIN
Changes since 2.49: +17 -10 lines
Log Message:
Extensive changes to enable rtrace -oTW option for tracking ray contributions

File Contents

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