ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/rt/source.c
Revision: 1.9
Committed: Thu Jun 8 09:35:43 1989 UTC (34 years, 10 months ago) by greg
Content type: text/plain
Branch: MAIN
Changes since 1.8: +30 -25 lines
Log Message:
Changed sorting array to reduce memory exchange time in qsort

File Contents

# Content
1 /* Copyright (c) 1986 Regents of the University of California */
2
3 #ifndef lint
4 static char SCCSid[] = "$SunId$ LBL";
5 #endif
6
7 /*
8 * source.c - routines dealing with illumination sources.
9 *
10 * 8/20/85
11 */
12
13 #include "ray.h"
14
15 #include "octree.h"
16
17 #include "source.h"
18
19 #include "otypes.h"
20
21 #include "cone.h"
22
23 #include "face.h"
24
25 #include "random.h"
26
27
28 extern double dstrsrc; /* source distribution amount */
29 extern double shadthresh; /* relative shadow threshold */
30
31 SRCREC *source = NULL; /* our list of sources */
32 int nsources = 0; /* the number of sources */
33
34
35 marksources() /* find and mark source objects */
36 {
37 register OBJREC *o, *m;
38 register int i;
39
40 for (i = 0; i < nobjects; i++) {
41
42 o = objptr(i);
43
44 if (o->omod == OVOID)
45 continue;
46
47 m = objptr(o->omod);
48
49 if (m->otype != MAT_LIGHT &&
50 m->otype != MAT_ILLUM &&
51 m->otype != MAT_GLOW &&
52 m->otype != MAT_SPOT)
53 continue;
54
55 if (m->oargs.nfargs != (m->otype == MAT_GLOW ? 4 :
56 m->otype == MAT_SPOT ? 7 : 3))
57 objerror(m, USER, "bad # arguments");
58
59 if (m->otype == MAT_GLOW &&
60 o->otype != OBJ_SOURCE &&
61 m->oargs.farg[3] <= FTINY)
62 continue; /* don't bother */
63
64 if (source == NULL)
65 source = (SRCREC *)malloc(sizeof(SRCREC));
66 else
67 source = (SRCREC *)realloc((char *)source,
68 (unsigned)(nsources+1)*sizeof(SRCREC));
69 if (source == NULL)
70 error(SYSTEM, "out of memory in marksources");
71
72 newsource(&source[nsources], o);
73
74 if (m->otype == MAT_GLOW) {
75 source[nsources].sflags |= SPROX;
76 source[nsources].sl.prox = m->oargs.farg[3];
77 if (o->otype == OBJ_SOURCE)
78 source[nsources].sflags |= SSKIP;
79 } else if (m->otype == MAT_SPOT) {
80 source[nsources].sflags |= SSPOT;
81 source[nsources].sl.s = makespot(m);
82 }
83 nsources++;
84 }
85 }
86
87
88 newsource(src, so) /* add a source to the array */
89 register SRCREC *src;
90 register OBJREC *so;
91 {
92 double cos(), tan(), sqrt();
93 double theta;
94 FACE *f;
95 CONE *co;
96 int j;
97 register int i;
98
99 src->sflags = 0;
100 src->nhits = 1; src->ntests = 2; /* start probability = 1/2 */
101 src->so = so;
102
103 switch (so->otype) {
104 case OBJ_SOURCE:
105 if (so->oargs.nfargs != 4)
106 objerror(so, USER, "bad arguments");
107 src->sflags |= SDISTANT;
108 VCOPY(src->sloc, so->oargs.farg);
109 if (normalize(src->sloc) == 0.0)
110 objerror(so, USER, "zero direction");
111 theta = PI/180.0/2.0 * so->oargs.farg[3];
112 if (theta <= FTINY)
113 objerror(so, USER, "zero size");
114 src->ss = theta >= PI/4 ? 1.0 : tan(theta);
115 src->ss2 = 2.0*PI * (1.0 - cos(theta));
116 break;
117 case OBJ_SPHERE:
118 VCOPY(src->sloc, so->oargs.farg);
119 src->ss = so->oargs.farg[3];
120 src->ss2 = PI * src->ss * src->ss;
121 break;
122 case OBJ_FACE:
123 /* get the face */
124 f = getface(so);
125 /* find the center */
126 for (j = 0; j < 3; j++) {
127 src->sloc[j] = 0.0;
128 for (i = 0; i < f->nv; i++)
129 src->sloc[j] += VERTEX(f,i)[j];
130 src->sloc[j] /= f->nv;
131 }
132 if (!inface(src->sloc, f))
133 objerror(so, USER, "cannot hit center");
134 src->ss = sqrt(f->area / PI);
135 src->ss2 = f->area;
136 break;
137 case OBJ_RING:
138 /* get the ring */
139 co = getcone(so, 0);
140 VCOPY(src->sloc, CO_P0(co));
141 if (CO_R0(co) > 0.0)
142 objerror(so, USER, "cannot hit center");
143 src->ss = CO_R1(co);
144 src->ss2 = PI * src->ss * src->ss;
145 break;
146 default:
147 objerror(so, USER, "illegal material");
148 }
149 }
150
151
152 SPOT *
153 makespot(m) /* make a spotlight */
154 register OBJREC *m;
155 {
156 extern double cos();
157 register SPOT *ns;
158
159 if ((ns = (SPOT *)malloc(sizeof(SPOT))) == NULL)
160 error(SYSTEM, "out of memory in makespot");
161 ns->siz = 2.0*PI * (1.0 - cos(PI/180.0/2.0 * m->oargs.farg[3]));
162 VCOPY(ns->aim, m->oargs.farg+4);
163 if ((ns->flen = normalize(ns->aim)) == 0.0)
164 objerror(m, USER, "zero focus vector");
165 return(ns);
166 }
167
168
169 double
170 srcray(sr, r, sn) /* send a ray to a source, return domega */
171 register RAY *sr; /* returned source ray */
172 RAY *r; /* ray which hit object */
173 register int sn; /* source number */
174 {
175 register double *norm = NULL; /* plane normal */
176 double ddot; /* (distance times) cosine */
177 FVECT vd;
178 double d;
179 register int i;
180
181 if (source[sn].sflags & SSKIP)
182 return(0.0); /* skip this source */
183
184 rayorigin(sr, r, SHADOW, 1.0); /* ignore limits */
185
186 sr->rsrc = sn; /* remember source */
187 /* get source direction */
188 if (source[sn].sflags & SDISTANT)
189 /* constant direction */
190 VCOPY(sr->rdir, source[sn].sloc);
191 else { /* compute direction */
192 for (i = 0; i < 3; i++)
193 sr->rdir[i] = source[sn].sloc[i] - sr->rorg[i];
194
195 if (source[sn].so->otype == OBJ_FACE)
196 norm = getface(source[sn].so)->norm;
197 else if (source[sn].so->otype == OBJ_RING)
198 norm = getcone(source[sn].so,0)->ad;
199
200 if (norm != NULL && (ddot = -DOT(sr->rdir, norm)) <= FTINY)
201 return(0.0); /* behind surface! */
202 }
203 if (dstrsrc > FTINY) {
204 /* distribute source direction */
205 for (i = 0; i < 3; i++)
206 vd[i] = dstrsrc * source[sn].ss * (1.0 - 2.0*frandom());
207
208 if (norm != NULL) { /* project offset */
209 d = DOT(vd, norm);
210 for (i = 0; i < 3; i++)
211 vd[i] -= d * norm[i];
212 }
213 for (i = 0; i < 3; i++) /* offset source direction */
214 sr->rdir[i] += vd[i];
215
216 } else if (source[sn].sflags & SDISTANT)
217 /* already normalized */
218 return(source[sn].ss2);
219
220 if ((d = normalize(sr->rdir)) == 0.0)
221 /* at source! */
222 return(0.0);
223
224 if (source[sn].sflags & SDISTANT)
225 /* domega constant */
226 return(source[sn].ss2);
227
228 else {
229 /* check proximity */
230 if (source[sn].sflags & SPROX &&
231 d > source[sn].sl.prox)
232 return(0.0);
233
234 if (norm != NULL)
235 ddot /= d;
236 else
237 ddot = 1.0;
238 /* check angle */
239 if (source[sn].sflags & SSPOT) {
240 if (source[sn].sl.s->siz < 2.0*PI *
241 (1.0 + DOT(source[sn].sl.s->aim,sr->rdir)))
242 return(0.0);
243 d += source[sn].sl.s->flen;
244 }
245 /* return domega */
246 return(ddot*source[sn].ss2/(d*d));
247 }
248 }
249
250
251 sourcehit(r) /* check to see if ray hit distant source */
252 register RAY *r;
253 {
254 int first, last;
255 register int i;
256
257 if (r->rsrc >= 0) { /* check only one if aimed */
258 first = last = r->rsrc;
259 } else { /* otherwise check all */
260 first = 0; last = nsources-1;
261 }
262 for (i = first; i <= last; i++)
263 if (source[i].sflags & SDISTANT)
264 /*
265 * Check to see if ray is within
266 * solid angle of source.
267 */
268 if (2.0*PI * (1.0 - DOT(source[i].sloc,r->rdir))
269 <= source[i].ss2) {
270 r->ro = source[i].so;
271 if (!(source[i].sflags & SSKIP))
272 break;
273 }
274
275 if (r->ro != NULL) {
276 for (i = 0; i < 3; i++)
277 r->ron[i] = -r->rdir[i];
278 r->rod = 1.0;
279 r->rofs = 1.0; setident4(r->rofx);
280 r->robs = 1.0; setident4(r->robx);
281 return(1);
282 }
283 return(0);
284 }
285
286
287 static int
288 cntcmp(sc1, sc2) /* contribution compare (descending) */
289 register CNTPTR *sc1, *sc2;
290 {
291 if (sc1->brt > sc2->brt)
292 return(-1);
293 if (sc1->brt < sc2->brt)
294 return(1);
295 return(0);
296 }
297
298
299 direct(r, f, p) /* add direct component */
300 RAY *r; /* ray that hit surface */
301 int (*f)(); /* direct component coefficient function */
302 char *p; /* data for f */
303 {
304 register int sn;
305 register CONTRIB *srccnt;
306 register CNTPTR *cntord;
307 int ncnts;
308 double ourthresh, prob, hwt, test2, hit2;
309 RAY sr;
310
311 srccnt = (CONTRIB *)malloc(nsources*sizeof(CONTRIB));
312 cntord = (CNTPTR *)malloc(nsources*sizeof(CNTPTR));
313 if (srccnt == NULL || cntord == NULL)
314 error(SYSTEM, "out of memory in direct");
315 /* modify threshold */
316 ourthresh = shadthresh / r->rweight;
317 /* potential contributions */
318 for (sn = 0; sn < nsources; sn++) {
319 cntord[sn].sno = sn;
320 cntord[sn].brt = 0.0;
321 /* get source ray */
322 if ((srccnt[sn].dom = srcray(&sr, r, sn)) == 0.0)
323 continue;
324 VCOPY(srccnt[sn].dir, sr.rdir);
325 /* compute coefficient */
326 (*f)(srccnt[sn].val, p, srccnt[sn].dir, srccnt[sn].dom);
327 cntord[sn].brt = bright(srccnt[sn].val);
328 if (cntord[sn].brt <= FTINY)
329 continue;
330 /* compute intersection */
331 if (!( source[sn].sflags & SDISTANT ?
332 sourcehit(&sr) :
333 (*ofun[source[sn].so->otype].funp)
334 (source[sn].so, &sr) ))
335 continue;
336 /* compute contribution */
337 rayshade(&sr, sr.ro->omod);
338 multcolor(srccnt[sn].val, sr.rcol);
339 cntord[sn].brt = bright(srccnt[sn].val);
340 }
341 /* sort contributions */
342 qsort(cntord, nsources, sizeof(CNTPTR), cntcmp);
343 hit2 = 0.0; test2 = FTINY;
344 /* find last */
345 sn = 0; ncnts = nsources;
346 while (sn < ncnts-1) {
347 register int m;
348 m = (sn + ncnts) >> 1;
349 if (cntord[m].brt > 0.0)
350 sn = m;
351 else
352 ncnts = m;
353 }
354 /* accumulate tail */
355 for (sn = ncnts-1; sn > 0; sn--)
356 cntord[sn-1].brt += cntord[sn].brt;
357 /* shadow testing */
358 for (sn = 0; sn < ncnts; sn++) {
359 /* tail below threshold? */
360 if (cntord[sn].brt < ourthresh*bright(r->rcol))
361 break;
362 /* get statistics */
363 hwt = (double)source[cntord[sn].sno].nhits /
364 (double)source[cntord[sn].sno].ntests;
365 test2 += hwt;
366 source[cntord[sn].sno].ntests++;
367 /* test for hit */
368 rayorigin(&sr, r, SHADOW, 1.0);
369 VCOPY(sr.rdir, srccnt[cntord[sn].sno].dir);
370 if (localhit(&sr, &thescene) &&
371 sr.ro != source[cntord[sn].sno].so) {
372 /* check for transmission */
373 if (sr.clipset != NULL && inset(sr.clipset,sr.ro->omod))
374 raytrans(&sr); /* object is clipped */
375 else
376 rayshade(&sr, sr.ro->omod);
377 if (bright(sr.rcol) <= FTINY)
378 continue; /* missed! */
379 (*f)(srccnt[cntord[sn].sno].val, p,
380 srccnt[cntord[sn].sno].dir,
381 srccnt[cntord[sn].sno].dom);
382 multcolor(srccnt[cntord[sn].sno].val, sr.rcol);
383 }
384 /* add contribution if hit */
385 addcolor(r->rcol, srccnt[cntord[sn].sno].val);
386 hit2 += hwt;
387 source[cntord[sn].sno].nhits++;
388 }
389 /* weighted hit rate */
390 hwt = hit2 / test2;
391 #ifdef DEBUG
392 fprintf(stderr, "%d tested, %d untested, %f hit rate\n",
393 sn, ncnts-sn, hwt);
394 #endif
395 /* add in untested sources */
396 for ( ; sn < ncnts; sn++) {
397 prob = hwt * (double)source[cntord[sn].sno].nhits /
398 (double)source[cntord[sn].sno].ntests;
399 scalecolor(srccnt[cntord[sn].sno].val, prob);
400 addcolor(r->rcol, srccnt[cntord[sn].sno].val);
401 }
402 free(srccnt);
403 free(cntord);
404 }
405
406
407 #define wrongsource(m, r) (m->otype!=MAT_ILLUM && \
408 r->rsrc>=0 && \
409 source[r->rsrc].so!=r->ro)
410
411 #define badambient(m, r) ((r->crtype&(AMBIENT|SHADOW))==AMBIENT && \
412 !(r->rtype&REFLECTED) && /* hack! */\
413 !(m->otype==MAT_GLOW&&r->rot>m->oargs.farg[3]))
414
415 #define passillum(m, r) (m->otype==MAT_ILLUM && \
416 !(r->rsrc>=0&&source[r->rsrc].so==r->ro))
417
418
419 m_light(m, r) /* ray hit a light source */
420 register OBJREC *m;
421 register RAY *r;
422 {
423 /* check for behind */
424 if (r->rod < 0.0)
425 return;
426 /* check for over-counting */
427 if (wrongsource(m, r) || badambient(m, r))
428 return;
429 /* check for passed illum */
430 if (passillum(m, r)) {
431
432 if (m->oargs.nsargs < 1 || !strcmp(m->oargs.sarg[0], VOIDID))
433 raytrans(r);
434 else
435 rayshade(r, modifier(m->oargs.sarg[0]));
436
437 /* otherwise treat as source */
438 } else {
439 /* get distribution pattern */
440 raytexture(r, m->omod);
441 /* get source color */
442 setcolor(r->rcol, m->oargs.farg[0],
443 m->oargs.farg[1],
444 m->oargs.farg[2]);
445 /* modify value */
446 multcolor(r->rcol, r->pcol);
447 }
448 }
449
450
451 o_source() {} /* intersection with a source is done elsewhere */