ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/rt/source.c
Revision: 2.29
Committed: Sat Feb 22 02:07:29 2003 UTC (21 years, 2 months ago) by greg
Content type: text/plain
Branch: MAIN
Changes since 2.28: +86 -9 lines
Log Message:
Changes and check-in for 3.5 release
Includes new source files and modifications not recorded for many years
See ray/doc/notes/ReleaseNotes for notes between 3.1 and 3.5 release

File Contents

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