ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/rt/virtuals.c
Revision: 2.9
Committed: Wed Apr 23 00:52:34 2003 UTC (21 years ago) by greg
Content type: text/plain
Branch: MAIN
Changes since 2.8: +1 -1 lines
Log Message:
Added (void *) cast to realloc calls

File Contents

# Content
1 #ifndef lint
2 static const char RCSid[] = "$Id$";
3 #endif
4 /*
5 * Routines for simulating virtual light sources
6 * Thus far, we only support planar mirrors.
7 *
8 * External symbols declared in source.h
9 */
10
11 #include "copyright.h"
12
13 #include "ray.h"
14
15 #include "otypes.h"
16
17 #include "source.h"
18
19 #include "random.h"
20
21 #define MINSAMPLES 16 /* minimum number of pretest samples */
22 #define STESTMAX 32 /* maximum seeks per sample */
23
24
25 static OBJECT *vobject; /* virtual source objects */
26 static int nvobjects = 0; /* number of virtual source objects */
27
28
29 void
30 markvirtuals() /* find and mark virtual sources */
31 {
32 register OBJREC *o;
33 register int i;
34 /* check number of direct relays */
35 if (directrelay <= 0)
36 return;
37 /* find virtual source objects */
38 for (i = 0; i < nobjects; i++) {
39 o = objptr(i);
40 if (!issurface(o->otype) || o->omod == OVOID)
41 continue;
42 if (!isvlight(vsmaterial(o)->otype))
43 continue;
44 if (sfun[o->otype].of == NULL ||
45 sfun[o->otype].of->getpleq == NULL) {
46 objerror(o,WARNING,"secondary sources not supported");
47 continue;
48 }
49 if (nvobjects == 0)
50 vobject = (OBJECT *)malloc(sizeof(OBJECT));
51 else
52 vobject = (OBJECT *)realloc((void *)vobject,
53 (unsigned)(nvobjects+1)*sizeof(OBJECT));
54 if (vobject == NULL)
55 error(SYSTEM, "out of memory in addvirtuals");
56 vobject[nvobjects++] = i;
57 }
58 if (nvobjects == 0)
59 return;
60 #ifdef DEBUG
61 fprintf(stderr, "found %d virtual source objects\n", nvobjects);
62 #endif
63 /* append virtual sources */
64 for (i = nsources; i-- > 0; )
65 addvirtuals(i, directrelay);
66 /* done with our object list */
67 free((void *)vobject);
68 nvobjects = 0;
69 }
70
71
72 void
73 addvirtuals(sn, nr) /* add virtuals associated with source */
74 int sn;
75 int nr;
76 {
77 register int i;
78 /* check relay limit first */
79 if (nr <= 0)
80 return;
81 if (source[sn].sflags & SSKIP)
82 return;
83 /* check each virtual object for projection */
84 for (i = 0; i < nvobjects; i++)
85 /* vproject() calls us recursively */
86 vproject(objptr(vobject[i]), sn, nr-1);
87 }
88
89
90 void
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 = lastmod(objndx(m), m->oargs.sarg[0])) == OVOID)
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 = -1.;
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 stestlim, ssn;
305 int nhit, nok;
306 register int i, n;
307 /* return if pretesting disabled */
308 if (vspretest <= 0)
309 return(f);
310 /* get surface normal */
311 getplaneq(onorm, o);
312 /* set number of rays to sample */
313 if (source[sn].sflags & SDISTANT) {
314 /* 32. == heuristic constant */
315 n = 32.*or2/(thescene.cusize*thescene.cusize)*vspretest + .5;
316 } else {
317 for (i = 0; i < 3; i++)
318 offsdir[i] = source[sn].sloc[i] - oc[i];
319 d = DOT(offsdir,offsdir);
320 if (d <= FTINY)
321 n = 2.*PI * vspretest + .5;
322 else
323 n = 2.*PI * (1.-sqrt(1./(1.+or2/d)))*vspretest + .5;
324 }
325 if (n < MINSAMPLES) n = MINSAMPLES;
326 #ifdef DEBUG
327 fprintf(stderr, "pretesting source %d in object %s with %d rays\n",
328 sn, o->oname, n);
329 #endif
330 /* sample */
331 or = sqrt(or2);
332 stestlim = n*STESTMAX;
333 ssn = 0;
334 nhit = nok = 0;
335 initsrcindex(&si);
336 while (n-- > 0) {
337 /* get sample point */
338 do {
339 if (ssn >= stestlim) {
340 #ifdef DEBUG
341 fprintf(stderr, "\ttoo hard to hit\n");
342 #endif
343 return(f); /* too small a target! */
344 }
345 multisamp(offsdir, 3, urand(sn*931+5827+ssn));
346 for (i = 0; i < 3; i++)
347 offsdir[i] = or*(1. - 2.*offsdir[i]);
348 ssn++;
349 d = 1. - DOT(offsdir, onorm);
350 for (i = 0; i < 3; i++) {
351 sr.rorg[i] = oc[i] + offsdir[i] + d*onorm[i];
352 sr.rdir[i] = -onorm[i];
353 }
354 sr.rmax = 0.0;
355 rayorigin(&sr, NULL, PRIMARY, 1.0);
356 } while (!(*ofun[o->otype].funp)(o, &sr));
357 /* check against source */
358 VCOPY(sr.rorg, sr.rop); /* starting from intersection */
359 samplendx++;
360 if (si.sp >= si.np-1 ||
361 !srcray(&sr, NULL, &si) || sr.rsrc != sn) {
362 si.sn = sn-1; /* reset index to our source */
363 si.np = 0;
364 if (!srcray(&sr, NULL, &si) || sr.rsrc != sn)
365 continue; /* can't get there from here */
366 }
367 sr.revf = srcvalue;
368 rayvalue(&sr); /* check sample validity */
369 if (bright(sr.rcol) <= FTINY)
370 continue;
371 nok++; /* got sample; check obstructions */
372 rayclear(&sr);
373 sr.revf = raytrace;
374 rayvalue(&sr);
375 if (bright(sr.rcol) > FTINY)
376 nhit++;
377 if (nhit > 0 && nhit < nok) {
378 #ifdef DEBUG
379 fprintf(stderr, "\tpartially occluded\n");
380 #endif
381 return(f); /* need to shadow test */
382 }
383 }
384 if (nhit == 0) {
385 #ifdef DEBUG
386 fprintf(stderr, "\t0%% hit rate\n");
387 #endif
388 return(f | SSKIP); /* 0% hit rate: totally occluded */
389 }
390 #ifdef DEBUG
391 fprintf(stderr, "\t100%% hit rate\n");
392 #endif
393 return(f & ~SFOLLOW); /* 100% hit rate: no occlusion */
394 }
395
396
397 #ifdef DEBUG
398 void
399 virtverb(sn, fp) /* print verbose description of virtual source */
400 register int sn;
401 FILE *fp;
402 {
403 register int i;
404
405 fprintf(fp, "%s virtual source %d in %s %s\n",
406 source[sn].sflags & SDISTANT ? "distant" : "local",
407 sn, ofun[source[sn].so->otype].funame,
408 source[sn].so->oname);
409 fprintf(fp, "\tat (%f,%f,%f)\n",
410 source[sn].sloc[0], source[sn].sloc[1], source[sn].sloc[2]);
411 fprintf(fp, "\tlinked to source %d (%s)\n",
412 source[sn].sa.sv.sn, source[source[sn].sa.sv.sn].so->oname);
413 if (source[sn].sflags & SFOLLOW)
414 fprintf(fp, "\talways followed\n");
415 else
416 fprintf(fp, "\tnever followed\n");
417 if (!(source[sn].sflags & SSPOT))
418 return;
419 fprintf(fp, "\twith spot aim (%f,%f,%f) and size %f\n",
420 source[sn].sl.s->aim[0], source[sn].sl.s->aim[1],
421 source[sn].sl.s->aim[2], source[sn].sl.s->siz);
422 }
423 #endif