ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/rt/virtuals.c
Revision: 1.31
Committed: Fri Nov 8 16:34:27 1991 UTC (32 years, 5 months ago) by greg
Content type: text/plain
Branch: MAIN
Changes since 1.30: +19 -2 lines
Log Message:
fixed bug that disallowed illums to pass to virtual sources

File Contents

# Content
1 /* Copyright (c) 1991 Regents of the University of California */
2
3 #ifndef lint
4 static char SCCSid[] = "$SunId$ LBL";
5 #endif
6
7 /*
8 * Routines for simulating virtual light sources
9 * Thus far, we only support planar mirrors.
10 */
11
12 #include "ray.h"
13
14 #include "octree.h"
15
16 #include "otypes.h"
17
18 #include "source.h"
19
20 #include "random.h"
21
22 #define MINSAMPLES 16 /* minimum number of pretest samples */
23 #define STESTMAX 32 /* maximum seeks per sample */
24
25
26 double getdisk();
27
28 static OBJECT *vobject; /* virtual source objects */
29 static int nvobjects = 0; /* number of virtual source objects */
30
31
32 markvirtuals() /* find and mark virtual sources */
33 {
34 register OBJREC *o;
35 register int i;
36 /* check number of direct relays */
37 if (directrelay <= 0)
38 return;
39 /* find virtual source objects */
40 for (i = 0; i < nobjects; i++) {
41 o = objptr(i);
42 if (!issurface(o->otype) || o->omod == OVOID)
43 continue;
44 if (!isvlight(vsmaterial(o)->otype))
45 continue;
46 if (sfun[o->otype].of == NULL ||
47 sfun[o->otype].of->getpleq == NULL) {
48 objerror(o,WARNING,"secondary sources not supported");
49 continue;
50 }
51 if (nvobjects == 0)
52 vobject = (OBJECT *)malloc(sizeof(OBJECT));
53 else
54 vobject = (OBJECT *)realloc((char *)vobject,
55 (unsigned)(nvobjects+1)*sizeof(OBJECT));
56 if (vobject == NULL)
57 error(SYSTEM, "out of memory in addvirtuals");
58 vobject[nvobjects++] = i;
59 }
60 if (nvobjects == 0)
61 return;
62 #ifdef DEBUG
63 fprintf(stderr, "found %d virtual source objects\n", nvobjects);
64 #endif
65 /* append virtual sources */
66 for (i = nsources; i-- > 0; )
67 addvirtuals(i, directrelay);
68 /* done with our object list */
69 free((char *)vobject);
70 nvobjects = 0;
71 }
72
73
74 addvirtuals(sn, nr) /* add virtuals associated with source */
75 int sn;
76 int nr;
77 {
78 register int i;
79 /* check relay limit first */
80 if (nr <= 0)
81 return;
82 if (source[sn].sflags & SSKIP)
83 return;
84 /* check each virtual object for projection */
85 for (i = 0; i < nvobjects; i++)
86 /* vproject() calls us recursively */
87 vproject(objptr(vobject[i]), sn, nr-1);
88 }
89
90
91 vproject(o, sn, n) /* create projected source(s) if they exist */
92 OBJREC *o;
93 int sn;
94 int n;
95 {
96 register int i;
97 register VSMATERIAL *vsmat;
98 MAT4 proj;
99 int ns;
100
101 if (o == source[sn].so) /* objects cannot project themselves */
102 return;
103 /* get virtual source material */
104 vsmat = sfun[vsmaterial(o)->otype].mf;
105 /* project virtual sources */
106 for (i = 0; i < vsmat->nproj; i++)
107 if ((*vsmat->vproj)(proj, o, &source[sn], i))
108 if ((ns = makevsrc(o, sn, proj)) >= 0) {
109 source[ns].sa.sv.pn = i;
110 #ifdef DEBUG
111 virtverb(ns, stderr);
112 #endif
113 addvirtuals(ns, n);
114 }
115 }
116
117
118 OBJREC *
119 vsmaterial(o) /* get virtual source material pointer */
120 OBJREC *o;
121 {
122 register int i;
123 register OBJREC *m;
124
125 i = o->omod;
126 m = objptr(i);
127 if (m->otype != MAT_ILLUM || m->oargs.nsargs < 1 ||
128 !strcmp(m->oargs.sarg[0], VOIDID) ||
129 (i = modifier(m->oargs.sarg[0])) < 0)
130 return(m); /* direct modifier */
131 return(objptr(i)); /* illum alternate */
132 }
133
134
135 int
136 makevsrc(op, sn, pm) /* make virtual source if reasonable */
137 OBJREC *op;
138 register int sn;
139 MAT4 pm;
140 {
141 FVECT nsloc, nsnorm, ocent, v;
142 double maxrad2, d;
143 int nsflags;
144 SPOT theirspot, ourspot;
145 register int i;
146
147 nsflags = source[sn].sflags | (SVIRTUAL|SSPOT|SFOLLOW);
148 /* get object center and max. radius */
149 maxrad2 = getdisk(ocent, op, sn);
150 if (maxrad2 <= FTINY) /* too small? */
151 return(-1);
152 /* get location and spot */
153 if (source[sn].sflags & SDISTANT) { /* distant source */
154 if (source[sn].sflags & SPROX)
155 return(-1); /* should never get here! */
156 multv3(nsloc, source[sn].sloc, pm);
157 normalize(nsloc);
158 VCOPY(ourspot.aim, ocent);
159 ourspot.siz = PI*maxrad2;
160 ourspot.flen = 0.;
161 if (source[sn].sflags & SSPOT) {
162 multp3(theirspot.aim, source[sn].sl.s->aim, pm);
163 /* adjust for source size */
164 d = sqrt(dist2(ourspot.aim, theirspot.aim));
165 d = sqrt(source[sn].sl.s->siz/PI) + d*source[sn].srad;
166 theirspot.siz = PI*d*d;
167 ourspot.flen = theirspot.flen = source[sn].sl.s->flen;
168 d = ourspot.siz;
169 if (!commonbeam(&ourspot, &theirspot, nsloc))
170 return(-1); /* no overlap */
171 if (ourspot.siz < d-FTINY) { /* it shrunk */
172 d = beamdisk(v, op, &ourspot, nsloc);
173 if (d <= FTINY)
174 return(-1);
175 if (d < maxrad2) {
176 maxrad2 = d;
177 VCOPY(ocent, v);
178 }
179 }
180 }
181 } else { /* local source */
182 multp3(nsloc, source[sn].sloc, pm);
183 for (i = 0; i < 3; i++)
184 ourspot.aim[i] = ocent[i] - nsloc[i];
185 if ((d = normalize(ourspot.aim)) == 0.)
186 return(-1); /* at source!! */
187 if (source[sn].sflags & SPROX && d > source[sn].sl.prox)
188 return(-1); /* too far away */
189 ourspot.flen = 0.;
190 /* adjust for source size */
191 d = (sqrt(maxrad2) + source[sn].srad) / d;
192 if (d < 1.-FTINY)
193 ourspot.siz = 2.*PI*(1. - sqrt(1.-d*d));
194 else
195 nsflags &= ~SSPOT;
196 if (source[sn].sflags & SSPOT) {
197 copystruct(&theirspot, source[sn].sl.s);
198 multv3(theirspot.aim, source[sn].sl.s->aim, pm);
199 normalize(theirspot.aim);
200 if (nsflags & SSPOT) {
201 ourspot.flen = theirspot.flen;
202 d = ourspot.siz;
203 if (!commonspot(&ourspot, &theirspot, nsloc))
204 return(-1); /* no overlap */
205 } else {
206 nsflags |= SSPOT;
207 copystruct(&ourspot, &theirspot);
208 d = 2.*ourspot.siz;
209 }
210 if (ourspot.siz < d-FTINY) { /* it shrunk */
211 d = spotdisk(v, op, &ourspot, nsloc);
212 if (d <= FTINY)
213 return(-1);
214 if (d < maxrad2) {
215 maxrad2 = d;
216 VCOPY(ocent, v);
217 }
218 }
219 }
220 if (source[sn].sflags & SFLAT) { /* behind source? */
221 multv3(nsnorm, source[sn].snorm, pm);
222 normalize(nsnorm);
223 if (nsflags & SSPOT && !checkspot(&ourspot, nsnorm))
224 return(-1);
225 }
226 }
227 /* pretest visibility */
228 nsflags = vstestvis(nsflags, op, ocent, maxrad2, sn);
229 if (nsflags & SSKIP)
230 return(-1); /* obstructed */
231 /* it all checks out, so make it */
232 if ((i = newsource()) < 0)
233 goto memerr;
234 source[i].sflags = nsflags;
235 VCOPY(source[i].sloc, nsloc);
236 multv3(source[i].ss[SU], source[sn].ss[SU], pm);
237 multv3(source[i].ss[SV], source[sn].ss[SV], pm);
238 if (nsflags & SFLAT)
239 VCOPY(source[i].snorm, nsnorm);
240 else
241 multv3(source[i].ss[SW], source[sn].ss[SW], pm);
242 source[i].srad = source[sn].srad;
243 source[i].ss2 = source[sn].ss2;
244 if (nsflags & SSPOT) {
245 if ((source[i].sl.s = (SPOT *)malloc(sizeof(SPOT))) == NULL)
246 goto memerr;
247 copystruct(source[i].sl.s, &ourspot);
248 }
249 if (nsflags & SPROX)
250 source[i].sl.prox = source[sn].sl.prox;
251 source[i].sa.sv.sn = sn;
252 source[i].so = op;
253 return(i);
254 memerr:
255 error(SYSTEM, "out of memory in makevsrc");
256 }
257
258
259 double
260 getdisk(oc, op, sn) /* get visible object disk */
261 FVECT oc;
262 OBJREC *op;
263 register int sn;
264 {
265 double rad2, roffs, offs, d, rd, rdoto;
266 FVECT rnrm, nrm;
267 /* first, use object getdisk function */
268 rad2 = getmaxdisk(oc, op);
269 if (!(source[sn].sflags & SVIRTUAL))
270 return(rad2); /* all done for normal source */
271 /* check for correct side of relay surface */
272 roffs = getplaneq(rnrm, source[sn].so);
273 rd = DOT(rnrm, source[sn].sloc); /* source projection */
274 if (!(source[sn].sflags & SDISTANT))
275 rd -= roffs;
276 d = DOT(rnrm, oc) - roffs; /* disk distance to relay plane */
277 if ((d > 0.) ^ (rd > 0.))
278 return(rad2); /* OK if opposite sides */
279 if (d*d >= rad2)
280 return(0.); /* no relay is possible */
281 /* we need a closer look */
282 offs = getplaneq(nrm, op);
283 rdoto = DOT(rnrm, nrm);
284 if (d*d >= rad2*(1.-rdoto*rdoto))
285 return(0.); /* disk entirely on projection side */
286 /* should shrink disk but I'm lazy */
287 return(rad2);
288 }
289
290
291 int
292 vstestvis(f, o, oc, or2, sn) /* pretest source visibility */
293 int f; /* virtual source flags */
294 OBJREC *o; /* relay object */
295 FVECT oc; /* relay object center */
296 double or2; /* relay object radius squared */
297 register int sn; /* target source number */
298 {
299 RAY sr;
300 FVECT onorm;
301 FVECT offsdir;
302 SRCINDEX si;
303 double or, d;
304 int infront;
305 int stestlim, ssn;
306 int nhit, nok;
307 register int i, n;
308 /* return if pretesting disabled */
309 if (vspretest <= 0)
310 return(f);
311 /* get surface normal */
312 getplaneq(onorm, o);
313 /* set number of rays to sample */
314 if (source[sn].sflags & SDISTANT) {
315 /* 32. == heuristic constant */
316 n = 32.*or2/(thescene.cusize*thescene.cusize)*vspretest + .5;
317 infront = DOT(onorm, source[sn].sloc) > 0.;
318 } else {
319 for (i = 0; i < 3; i++)
320 offsdir[i] = source[sn].sloc[i] - oc[i];
321 d = DOT(offsdir,offsdir);
322 if (d <= FTINY)
323 n = 2.*PI * vspretest + .5;
324 else
325 n = 2.*PI * (1.-sqrt(1./(1.+or2/d)))*vspretest + .5;
326 infront = DOT(onorm, offsdir) > 0.;
327 }
328 if (n < MINSAMPLES) n = MINSAMPLES;
329 #ifdef DEBUG
330 fprintf(stderr, "pretesting source %d in object %s with %d rays\n",
331 sn, o->oname, n);
332 #endif
333 /* sample */
334 or = sqrt(or2);
335 stestlim = n*STESTMAX;
336 ssn = 0;
337 nhit = nok = 0;
338 while (n-- > 0) {
339 /* get sample point */
340 do {
341 if (ssn >= stestlim) {
342 #ifdef DEBUG
343 fprintf(stderr, "\ttoo hard to hit\n");
344 #endif
345 return(f); /* too small a target! */
346 }
347 multisamp(offsdir, 3, urand(sn*931+5827+ssn));
348 for (i = 0; i < 3; i++)
349 offsdir[i] = or*(1. - 2.*offsdir[i]);
350 ssn++;
351 for (i = 0; i < 3; i++)
352 sr.rorg[i] = oc[i] + offsdir[i];
353 d = DOT(offsdir,onorm);
354 if (infront)
355 for (i = 0; i < 3; i++) {
356 sr.rorg[i] -= (d-.0001)*onorm[i];
357 sr.rdir[i] = -onorm[i];
358 }
359 else
360 for (i = 0; i < 3; i++) {
361 sr.rorg[i] -= (d+.0001)*onorm[i];
362 sr.rdir[i] = onorm[i];
363 }
364 rayorigin(&sr, NULL, PRIMARY, 1.0);
365 } while (!(*ofun[o->otype].funp)(o, &sr));
366 /* check against source */
367 initsrcindex(&si);
368 si.sn = sn;
369 nopart(&si, &sr);
370 samplendx++;
371 if (!srcray(&sr, NULL, &si) || sr.rsrc != sn)
372 continue;
373 sr.revf = srcvalue;
374 rayvalue(&sr);
375 if (bright(sr.rcol) <= FTINY)
376 continue;
377 nok++;
378 /* check against obstructions */
379 rayclear(&sr);
380 sr.revf = raytrace;
381 rayvalue(&sr);
382 if (bright(sr.rcol) > FTINY)
383 nhit++;
384 if (nhit > 0 && nhit < nok) {
385 #ifdef DEBUG
386 fprintf(stderr, "\tpartially occluded\n");
387 #endif
388 return(f); /* need to shadow test */
389 }
390 }
391 if (nhit == 0) {
392 #ifdef DEBUG
393 fprintf(stderr, "\t0%% hit rate\n");
394 #endif
395 return(f | SSKIP); /* 0% hit rate: totally occluded */
396 }
397 #ifdef DEBUG
398 fprintf(stderr, "\t100%% hit rate\n");
399 #endif
400 return(f & ~SFOLLOW); /* 100% hit rate: no occlusion */
401 }
402
403
404 #ifdef DEBUG
405 virtverb(sn, fp) /* print verbose description of virtual source */
406 register int sn;
407 FILE *fp;
408 {
409 register int i;
410
411 fprintf(fp, "%s virtual source %d in %s %s\n",
412 source[sn].sflags & SDISTANT ? "distant" : "local",
413 sn, ofun[source[sn].so->otype].funame,
414 source[sn].so->oname);
415 fprintf(fp, "\tat (%f,%f,%f)\n",
416 source[sn].sloc[0], source[sn].sloc[1], source[sn].sloc[2]);
417 fprintf(fp, "\tlinked to source %d (%s)\n",
418 source[sn].sa.sv.sn, source[source[sn].sa.sv.sn].so->oname);
419 if (source[sn].sflags & SFOLLOW)
420 fprintf(fp, "\talways followed\n");
421 else
422 fprintf(fp, "\tnever followed\n");
423 if (!(source[sn].sflags & SSPOT))
424 return;
425 fprintf(fp, "\twith spot aim (%f,%f,%f) and size %f\n",
426 source[sn].sl.s->aim[0], source[sn].sl.s->aim[1],
427 source[sn].sl.s->aim[2], source[sn].sl.s->siz);
428 }
429 #endif