ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/rt/virtuals.c
Revision: 2.7
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.6: +65 -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

# User Rev Content
1 greg 1.1 #ifndef lint
2 greg 2.7 static const char RCSid[] = "$Id$";
3 greg 1.1 #endif
4     /*
5     * Routines for simulating virtual light sources
6     * Thus far, we only support planar mirrors.
7 greg 2.7 *
8     * External symbols declared in source.h
9     */
10    
11     /* ====================================================================
12     * The Radiance Software License, Version 1.0
13     *
14     * Copyright (c) 1990 - 2002 The Regents of the University of California,
15     * through Lawrence Berkeley National Laboratory. All rights reserved.
16     *
17     * Redistribution and use in source and binary forms, with or without
18     * modification, are permitted provided that the following conditions
19     * are met:
20     *
21     * 1. Redistributions of source code must retain the above copyright
22     * notice, this list of conditions and the following disclaimer.
23     *
24     * 2. Redistributions in binary form must reproduce the above copyright
25     * notice, this list of conditions and the following disclaimer in
26     * the documentation and/or other materials provided with the
27     * distribution.
28     *
29     * 3. The end-user documentation included with the redistribution,
30     * if any, must include the following acknowledgment:
31     * "This product includes Radiance software
32     * (http://radsite.lbl.gov/)
33     * developed by the Lawrence Berkeley National Laboratory
34     * (http://www.lbl.gov/)."
35     * Alternately, this acknowledgment may appear in the software itself,
36     * if and wherever such third-party acknowledgments normally appear.
37     *
38     * 4. The names "Radiance," "Lawrence Berkeley National Laboratory"
39     * and "The Regents of the University of California" must
40     * not be used to endorse or promote products derived from this
41     * software without prior written permission. For written
42     * permission, please contact [email protected].
43     *
44     * 5. Products derived from this software may not be called "Radiance",
45     * nor may "Radiance" appear in their name, without prior written
46     * permission of Lawrence Berkeley National Laboratory.
47     *
48     * THIS SOFTWARE IS PROVIDED ``AS IS'' AND ANY EXPRESSED OR IMPLIED
49     * WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES
50     * OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
51     * DISCLAIMED. IN NO EVENT SHALL Lawrence Berkeley National Laboratory OR
52     * ITS CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
53     * SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
54     * LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF
55     * USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
56     * ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
57     * OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT
58     * OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
59     * SUCH DAMAGE.
60     * ====================================================================
61     *
62     * This software consists of voluntary contributions made by many
63     * individuals on behalf of Lawrence Berkeley National Laboratory. For more
64     * information on Lawrence Berkeley National Laboratory, please see
65     * <http://www.lbl.gov/>.
66 greg 1.1 */
67    
68     #include "ray.h"
69    
70 greg 1.3 #include "otypes.h"
71    
72 greg 1.1 #include "source.h"
73    
74 greg 1.7 #include "random.h"
75 greg 1.1
76 greg 1.22 #define MINSAMPLES 16 /* minimum number of pretest samples */
77     #define STESTMAX 32 /* maximum seeks per sample */
78 greg 1.1
79 greg 1.13
80 greg 1.1 static OBJECT *vobject; /* virtual source objects */
81     static int nvobjects = 0; /* number of virtual source objects */
82    
83    
84 greg 2.7 void
85 greg 1.1 markvirtuals() /* find and mark virtual sources */
86     {
87     register OBJREC *o;
88     register int i;
89     /* check number of direct relays */
90     if (directrelay <= 0)
91     return;
92     /* find virtual source objects */
93     for (i = 0; i < nobjects; i++) {
94     o = objptr(i);
95 greg 1.3 if (!issurface(o->otype) || o->omod == OVOID)
96 greg 1.1 continue;
97 greg 1.31 if (!isvlight(vsmaterial(o)->otype))
98 greg 1.1 continue;
99 greg 1.3 if (sfun[o->otype].of == NULL ||
100 greg 1.21 sfun[o->otype].of->getpleq == NULL) {
101     objerror(o,WARNING,"secondary sources not supported");
102     continue;
103     }
104 greg 1.1 if (nvobjects == 0)
105     vobject = (OBJECT *)malloc(sizeof(OBJECT));
106     else
107     vobject = (OBJECT *)realloc((char *)vobject,
108     (unsigned)(nvobjects+1)*sizeof(OBJECT));
109     if (vobject == NULL)
110     error(SYSTEM, "out of memory in addvirtuals");
111     vobject[nvobjects++] = i;
112     }
113     if (nvobjects == 0)
114     return;
115 greg 1.4 #ifdef DEBUG
116     fprintf(stderr, "found %d virtual source objects\n", nvobjects);
117     #endif
118 greg 1.1 /* append virtual sources */
119     for (i = nsources; i-- > 0; )
120 greg 1.7 addvirtuals(i, directrelay);
121 greg 1.1 /* done with our object list */
122 greg 2.7 free((void *)vobject);
123 greg 1.1 nvobjects = 0;
124     }
125    
126    
127 greg 2.7 void
128 greg 1.4 addvirtuals(sn, nr) /* add virtuals associated with source */
129     int sn;
130 greg 1.1 int nr;
131     {
132     register int i;
133     /* check relay limit first */
134     if (nr <= 0)
135     return;
136 greg 1.7 if (source[sn].sflags & SSKIP)
137     return;
138 greg 1.1 /* check each virtual object for projection */
139     for (i = 0; i < nvobjects; i++)
140 greg 1.3 /* vproject() calls us recursively */
141 greg 1.4 vproject(objptr(vobject[i]), sn, nr-1);
142 greg 1.1 }
143    
144    
145 greg 2.7 void
146 greg 1.4 vproject(o, sn, n) /* create projected source(s) if they exist */
147 greg 1.3 OBJREC *o;
148 greg 1.4 int sn;
149 greg 1.3 int n;
150     {
151     register int i;
152     register VSMATERIAL *vsmat;
153     MAT4 proj;
154 greg 1.4 int ns;
155    
156     if (o == source[sn].so) /* objects cannot project themselves */
157     return;
158 greg 1.3 /* get virtual source material */
159 greg 1.31 vsmat = sfun[vsmaterial(o)->otype].mf;
160 greg 1.3 /* project virtual sources */
161     for (i = 0; i < vsmat->nproj; i++)
162 greg 1.4 if ((*vsmat->vproj)(proj, o, &source[sn], i))
163     if ((ns = makevsrc(o, sn, proj)) >= 0) {
164 greg 1.17 source[ns].sa.sv.pn = i;
165 greg 1.4 #ifdef DEBUG
166 greg 1.6 virtverb(ns, stderr);
167 greg 1.4 #endif
168 greg 1.3 addvirtuals(ns, n);
169 greg 1.4 }
170 greg 1.31 }
171    
172    
173     OBJREC *
174     vsmaterial(o) /* get virtual source material pointer */
175     OBJREC *o;
176     {
177     register int i;
178     register OBJREC *m;
179    
180     i = o->omod;
181     m = objptr(i);
182     if (m->otype != MAT_ILLUM || m->oargs.nsargs < 1 ||
183     !strcmp(m->oargs.sarg[0], VOIDID) ||
184 gwlarson 2.6 (i = lastmod(objndx(m), m->oargs.sarg[0])) == OVOID)
185 greg 1.31 return(m); /* direct modifier */
186     return(objptr(i)); /* illum alternate */
187 greg 1.3 }
188    
189    
190 greg 1.4 int
191     makevsrc(op, sn, pm) /* make virtual source if reasonable */
192 greg 1.1 OBJREC *op;
193 greg 1.4 register int sn;
194 greg 1.1 MAT4 pm;
195     {
196 greg 1.9 FVECT nsloc, nsnorm, ocent, v;
197     double maxrad2, d;
198 greg 1.3 int nsflags;
199 greg 1.1 SPOT theirspot, ourspot;
200     register int i;
201 greg 1.3
202 greg 1.6 nsflags = source[sn].sflags | (SVIRTUAL|SSPOT|SFOLLOW);
203 greg 1.1 /* get object center and max. radius */
204 greg 1.6 maxrad2 = getdisk(ocent, op, sn);
205     if (maxrad2 <= FTINY) /* too small? */
206     return(-1);
207 greg 1.1 /* get location and spot */
208 greg 1.4 if (source[sn].sflags & SDISTANT) { /* distant source */
209     if (source[sn].sflags & SPROX)
210 greg 1.5 return(-1); /* should never get here! */
211 greg 1.4 multv3(nsloc, source[sn].sloc, pm);
212 greg 1.17 normalize(nsloc);
213 greg 1.6 VCOPY(ourspot.aim, ocent);
214     ourspot.siz = PI*maxrad2;
215 greg 2.5 ourspot.flen = -1.;
216 greg 1.4 if (source[sn].sflags & SSPOT) {
217     multp3(theirspot.aim, source[sn].sl.s->aim, pm);
218 greg 1.29 /* adjust for source size */
219 greg 1.19 d = sqrt(dist2(ourspot.aim, theirspot.aim));
220 greg 1.28 d = sqrt(source[sn].sl.s->siz/PI) + d*source[sn].srad;
221 greg 1.19 theirspot.siz = PI*d*d;
222     ourspot.flen = theirspot.flen = source[sn].sl.s->flen;
223 greg 1.9 d = ourspot.siz;
224 greg 1.6 if (!commonbeam(&ourspot, &theirspot, nsloc))
225 greg 1.9 return(-1); /* no overlap */
226     if (ourspot.siz < d-FTINY) { /* it shrunk */
227     d = beamdisk(v, op, &ourspot, nsloc);
228     if (d <= FTINY)
229     return(-1);
230     if (d < maxrad2) {
231     maxrad2 = d;
232     VCOPY(ocent, v);
233     }
234     }
235 greg 1.1 }
236     } else { /* local source */
237 greg 1.4 multp3(nsloc, source[sn].sloc, pm);
238 greg 1.6 for (i = 0; i < 3; i++)
239     ourspot.aim[i] = ocent[i] - nsloc[i];
240 greg 1.9 if ((d = normalize(ourspot.aim)) == 0.)
241 greg 1.6 return(-1); /* at source!! */
242 greg 1.9 if (source[sn].sflags & SPROX && d > source[sn].sl.prox)
243 greg 1.6 return(-1); /* too far away */
244     ourspot.flen = 0.;
245 greg 1.29 /* adjust for source size */
246 greg 1.28 d = (sqrt(maxrad2) + source[sn].srad) / d;
247 greg 1.19 if (d < 1.-FTINY)
248     ourspot.siz = 2.*PI*(1. - sqrt(1.-d*d));
249 greg 1.14 else
250     nsflags &= ~SSPOT;
251 greg 1.4 if (source[sn].sflags & SSPOT) {
252     copystruct(&theirspot, source[sn].sl.s);
253     multv3(theirspot.aim, source[sn].sl.s->aim, pm);
254 greg 1.17 normalize(theirspot.aim);
255 greg 1.14 if (nsflags & SSPOT) {
256     ourspot.flen = theirspot.flen;
257     d = ourspot.siz;
258     if (!commonspot(&ourspot, &theirspot, nsloc))
259     return(-1); /* no overlap */
260     } else {
261     nsflags |= SSPOT;
262     copystruct(&ourspot, &theirspot);
263     d = 2.*ourspot.siz;
264     }
265 greg 1.9 if (ourspot.siz < d-FTINY) { /* it shrunk */
266     d = spotdisk(v, op, &ourspot, nsloc);
267     if (d <= FTINY)
268     return(-1);
269     if (d < maxrad2) {
270     maxrad2 = d;
271     VCOPY(ocent, v);
272     }
273     }
274 greg 1.1 }
275 greg 1.4 if (source[sn].sflags & SFLAT) { /* behind source? */
276     multv3(nsnorm, source[sn].snorm, pm);
277 greg 1.17 normalize(nsnorm);
278 greg 1.20 if (nsflags & SSPOT && !checkspot(&ourspot, nsnorm))
279 greg 1.5 return(-1);
280 greg 1.1 }
281     }
282 greg 1.7 /* pretest visibility */
283     nsflags = vstestvis(nsflags, op, ocent, maxrad2, sn);
284     if (nsflags & SSKIP)
285     return(-1); /* obstructed */
286     /* it all checks out, so make it */
287 greg 1.6 if ((i = newsource()) < 0)
288 greg 1.1 goto memerr;
289 greg 1.6 source[i].sflags = nsflags;
290     VCOPY(source[i].sloc, nsloc);
291 greg 1.28 multv3(source[i].ss[SU], source[sn].ss[SU], pm);
292     multv3(source[i].ss[SV], source[sn].ss[SV], pm);
293 greg 1.3 if (nsflags & SFLAT)
294 greg 1.6 VCOPY(source[i].snorm, nsnorm);
295 greg 1.28 else
296     multv3(source[i].ss[SW], source[sn].ss[SW], pm);
297 greg 1.29 source[i].srad = source[sn].srad;
298 greg 1.28 source[i].ss2 = source[sn].ss2;
299 greg 1.14 if (nsflags & SSPOT) {
300     if ((source[i].sl.s = (SPOT *)malloc(sizeof(SPOT))) == NULL)
301     goto memerr;
302     copystruct(source[i].sl.s, &ourspot);
303     }
304 greg 1.3 if (nsflags & SPROX)
305 greg 1.6 source[i].sl.prox = source[sn].sl.prox;
306 greg 1.17 source[i].sa.sv.sn = sn;
307 greg 1.6 source[i].so = op;
308     return(i);
309 greg 1.1 memerr:
310     error(SYSTEM, "out of memory in makevsrc");
311     }
312    
313    
314 greg 1.6 double
315     getdisk(oc, op, sn) /* get visible object disk */
316     FVECT oc;
317     OBJREC *op;
318     register int sn;
319     {
320     double rad2, roffs, offs, d, rd, rdoto;
321     FVECT rnrm, nrm;
322     /* first, use object getdisk function */
323 greg 1.9 rad2 = getmaxdisk(oc, op);
324 greg 1.6 if (!(source[sn].sflags & SVIRTUAL))
325     return(rad2); /* all done for normal source */
326     /* check for correct side of relay surface */
327 greg 1.9 roffs = getplaneq(rnrm, source[sn].so);
328 greg 1.6 rd = DOT(rnrm, source[sn].sloc); /* source projection */
329     if (!(source[sn].sflags & SDISTANT))
330     rd -= roffs;
331     d = DOT(rnrm, oc) - roffs; /* disk distance to relay plane */
332     if ((d > 0.) ^ (rd > 0.))
333     return(rad2); /* OK if opposite sides */
334     if (d*d >= rad2)
335 greg 1.9 return(0.); /* no relay is possible */
336 greg 1.6 /* we need a closer look */
337 greg 1.9 offs = getplaneq(nrm, op);
338 greg 1.6 rdoto = DOT(rnrm, nrm);
339     if (d*d >= rad2*(1.-rdoto*rdoto))
340     return(0.); /* disk entirely on projection side */
341     /* should shrink disk but I'm lazy */
342     return(rad2);
343     }
344    
345    
346 greg 1.7 int
347     vstestvis(f, o, oc, or2, sn) /* pretest source visibility */
348     int f; /* virtual source flags */
349     OBJREC *o; /* relay object */
350     FVECT oc; /* relay object center */
351     double or2; /* relay object radius squared */
352     register int sn; /* target source number */
353 greg 1.1 {
354 greg 1.7 RAY sr;
355     FVECT onorm;
356     FVECT offsdir;
357 greg 1.28 SRCINDEX si;
358 greg 1.7 double or, d;
359 greg 1.16 int stestlim, ssn;
360 greg 1.11 int nhit, nok;
361 greg 1.7 register int i, n;
362     /* return if pretesting disabled */
363     if (vspretest <= 0)
364     return(f);
365     /* get surface normal */
366 greg 1.9 getplaneq(onorm, o);
367 greg 1.7 /* set number of rays to sample */
368 greg 1.8 if (source[sn].sflags & SDISTANT) {
369 greg 1.26 /* 32. == heuristic constant */
370     n = 32.*or2/(thescene.cusize*thescene.cusize)*vspretest + .5;
371 greg 1.8 } else {
372     for (i = 0; i < 3; i++)
373     offsdir[i] = source[sn].sloc[i] - oc[i];
374 greg 1.20 d = DOT(offsdir,offsdir);
375     if (d <= FTINY)
376     n = 2.*PI * vspretest + .5;
377     else
378     n = 2.*PI * (1.-sqrt(1./(1.+or2/d)))*vspretest + .5;
379 greg 1.8 }
380 greg 1.13 if (n < MINSAMPLES) n = MINSAMPLES;
381 greg 1.9 #ifdef DEBUG
382     fprintf(stderr, "pretesting source %d in object %s with %d rays\n",
383     sn, o->oname, n);
384     #endif
385 greg 1.7 /* sample */
386 greg 1.8 or = sqrt(or2);
387 greg 1.16 stestlim = n*STESTMAX;
388     ssn = 0;
389 greg 1.11 nhit = nok = 0;
390 greg 2.4 initsrcindex(&si);
391 greg 1.7 while (n-- > 0) {
392 greg 1.8 /* get sample point */
393     do {
394 greg 1.16 if (ssn >= stestlim) {
395 greg 1.9 #ifdef DEBUG
396     fprintf(stderr, "\ttoo hard to hit\n");
397     #endif
398 greg 1.8 return(f); /* too small a target! */
399 greg 1.9 }
400 greg 1.25 multisamp(offsdir, 3, urand(sn*931+5827+ssn));
401 greg 1.8 for (i = 0; i < 3; i++)
402 greg 1.23 offsdir[i] = or*(1. - 2.*offsdir[i]);
403 greg 1.16 ssn++;
404 greg 2.4 d = 1. - DOT(offsdir, onorm);
405     for (i = 0; i < 3; i++) {
406     sr.rorg[i] = oc[i] + offsdir[i] + d*onorm[i];
407     sr.rdir[i] = -onorm[i];
408     }
409 greg 2.3 sr.rmax = 0.0;
410 greg 1.8 rayorigin(&sr, NULL, PRIMARY, 1.0);
411     } while (!(*ofun[o->otype].funp)(o, &sr));
412     /* check against source */
413 greg 2.4 VCOPY(sr.rorg, sr.rop); /* starting from intersection */
414 greg 1.7 samplendx++;
415 greg 2.4 if (si.sp >= si.np-1 ||
416     !srcray(&sr, NULL, &si) || sr.rsrc != sn) {
417     si.sn = sn-1; /* reset index to our source */
418     si.np = 0;
419     if (!srcray(&sr, NULL, &si) || sr.rsrc != sn)
420     continue; /* can't get there from here */
421     }
422 greg 1.7 sr.revf = srcvalue;
423 greg 2.4 rayvalue(&sr); /* check sample validity */
424 greg 1.7 if (bright(sr.rcol) <= FTINY)
425     continue;
426 greg 2.4 nok++; /* got sample; check obstructions */
427 greg 1.18 rayclear(&sr);
428 greg 1.27 sr.revf = raytrace;
429     rayvalue(&sr);
430     if (bright(sr.rcol) > FTINY)
431 greg 1.11 nhit++;
432     if (nhit > 0 && nhit < nok) {
433 greg 1.9 #ifdef DEBUG
434 greg 1.11 fprintf(stderr, "\tpartially occluded\n");
435 greg 1.9 #endif
436 greg 1.11 return(f); /* need to shadow test */
437     }
438 greg 1.1 }
439 greg 1.9 if (nhit == 0) {
440     #ifdef DEBUG
441     fprintf(stderr, "\t0%% hit rate\n");
442     #endif
443 greg 1.7 return(f | SSKIP); /* 0% hit rate: totally occluded */
444 greg 1.9 }
445     #ifdef DEBUG
446     fprintf(stderr, "\t100%% hit rate\n");
447     #endif
448     return(f & ~SFOLLOW); /* 100% hit rate: no occlusion */
449 greg 1.1 }
450 greg 1.7
451 greg 1.4
452     #ifdef DEBUG
453 greg 2.7 void
454 greg 1.6 virtverb(sn, fp) /* print verbose description of virtual source */
455     register int sn;
456 greg 1.4 FILE *fp;
457     {
458     register int i;
459    
460     fprintf(fp, "%s virtual source %d in %s %s\n",
461 greg 1.6 source[sn].sflags & SDISTANT ? "distant" : "local",
462     sn, ofun[source[sn].so->otype].funame,
463     source[sn].so->oname);
464 greg 1.4 fprintf(fp, "\tat (%f,%f,%f)\n",
465 greg 1.6 source[sn].sloc[0], source[sn].sloc[1], source[sn].sloc[2]);
466 greg 1.4 fprintf(fp, "\tlinked to source %d (%s)\n",
467 greg 1.17 source[sn].sa.sv.sn, source[source[sn].sa.sv.sn].so->oname);
468 greg 1.6 if (source[sn].sflags & SFOLLOW)
469 greg 1.4 fprintf(fp, "\talways followed\n");
470     else
471     fprintf(fp, "\tnever followed\n");
472 greg 1.6 if (!(source[sn].sflags & SSPOT))
473 greg 1.4 return;
474     fprintf(fp, "\twith spot aim (%f,%f,%f) and size %f\n",
475 greg 1.6 source[sn].sl.s->aim[0], source[sn].sl.s->aim[1],
476     source[sn].sl.s->aim[2], source[sn].sl.s->siz);
477 greg 1.4 }
478     #endif