ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/rt/m_bsdf.c
Revision: 2.32
Committed: Sat Feb 18 19:12:49 2017 UTC (7 years, 2 months ago) by greg
Content type: text/plain
Branch: MAIN
Changes since 2.31: +44 -8 lines
Log Message:
Corrected double-counting of diffuse contributions

File Contents

# Content
1 #ifndef lint
2 static const char RCSid[] = "$Id: m_bsdf.c,v 2.31 2017/02/17 23:24:56 greg Exp $";
3 #endif
4 /*
5 * Shading for materials with BSDFs taken from XML data files
6 */
7
8 #include "copyright.h"
9
10 #include "ray.h"
11 #include "ambient.h"
12 #include "source.h"
13 #include "func.h"
14 #include "bsdf.h"
15 #include "random.h"
16 #include "pmapmat.h"
17
18 /*
19 * Arguments to this material include optional diffuse colors.
20 * String arguments include the BSDF and function files.
21 * A non-zero thickness causes the strange but useful behavior
22 * of translating transmitted rays this distance beneath the surface
23 * (opposite the surface normal) to bypass any intervening geometry.
24 * Translation only affects scattered, non-source-directed samples.
25 * A non-zero thickness has the further side-effect that an unscattered
26 * (view) ray will pass right through our material if it has any
27 * non-diffuse transmission, making the BSDF surface invisible. This
28 * shows the proxied geometry instead. Thickness has the further
29 * effect of turning off reflection on the hidden side so that rays
30 * heading in the opposite direction pass unimpeded through the BSDF
31 * surface. A paired surface may be placed on the opposide side of
32 * the detail geometry, less than this thickness away, if a two-way
33 * proxy is desired. Note that the sign of the thickness is important.
34 * A positive thickness hides geometry behind the BSDF surface and uses
35 * front reflectance and transmission properties. A negative thickness
36 * hides geometry in front of the surface when rays hit from behind,
37 * and applies only the transmission and backside reflectance properties.
38 * Reflection is ignored on the hidden side, as those rays pass through.
39 * The "up" vector for the BSDF is given by three variables, defined
40 * (along with the thickness) by the named function file, or '.' if none.
41 * Together with the surface normal, this defines the local coordinate
42 * system for the BSDF.
43 * We do not reorient the surface, so if the BSDF has no back-side
44 * reflectance and none is given in the real arguments, a BSDF surface
45 * with zero thickness will appear black when viewed from behind
46 * unless backface visibility is off.
47 * The diffuse arguments are added to components in the BSDF file,
48 * not multiplied. However, patterns affect this material as a multiplier
49 * on everything except non-diffuse reflection.
50 *
51 * Arguments for MAT_BSDF are:
52 * 6+ thick BSDFfile ux uy uz funcfile transform
53 * 0
54 * 0|3|6|9 rdf gdf bdf
55 * rdb gdb bdb
56 * rdt gdt bdt
57 */
58
59 /*
60 * Note that our reverse ray-tracing process means that the positions
61 * of incoming and outgoing vectors may be reversed in our calls
62 * to the BSDF library. This is fine, since the bidirectional nature
63 * of the BSDF (that's what the 'B' stands for) means it all works out.
64 */
65
66 typedef struct {
67 OBJREC *mp; /* material pointer */
68 RAY *pr; /* intersected ray */
69 FVECT pnorm; /* perturbed surface normal */
70 FVECT vray; /* local outgoing (return) vector */
71 double sr_vpsa[2]; /* sqrt of BSDF projected solid angle extrema */
72 RREAL toloc[3][3]; /* world to local BSDF coords */
73 RREAL fromloc[3][3]; /* local BSDF coords to world */
74 double thick; /* surface thickness */
75 SDData *sd; /* loaded BSDF data */
76 COLOR rdiff; /* diffuse reflection */
77 COLOR tdiff; /* diffuse transmission */
78 } BSDFDAT; /* BSDF material data */
79
80 #define cvt_sdcolor(cv, svp) ccy2rgb(&(svp)->spec, (svp)->cieY, cv)
81
82 /* Jitter ray sample according to projected solid angle and specjitter */
83 static void
84 bsdf_jitter(FVECT vres, BSDFDAT *ndp, double sr_psa)
85 {
86 VCOPY(vres, ndp->vray);
87 if (specjitter < 1.)
88 sr_psa *= specjitter;
89 if (sr_psa <= FTINY)
90 return;
91 vres[0] += sr_psa*(.5 - frandom());
92 vres[1] += sr_psa*(.5 - frandom());
93 normalize(vres);
94 }
95
96 /* Evaluate BSDF for direct component, returning true if OK to proceed */
97 static int
98 direct_bsdf_OK(COLOR cval, FVECT ldir, double omega, BSDFDAT *ndp)
99 {
100 int nsamp, ok = 0;
101 FVECT vsrc, vsmp, vjit;
102 double tomega;
103 double sf, tsr, sd[2];
104 COLOR csmp, cdiff;
105 double diffY;
106 SDValue sv;
107 SDError ec;
108 int i;
109 /* transform source direction */
110 if (SDmapDir(vsrc, ndp->toloc, ldir) != SDEnone)
111 return(0);
112 /* will discount diffuse portion */
113 switch ((vsrc[2] > 0)<<1 | (ndp->vray[2] > 0)) {
114 case 3:
115 if (ndp->sd->rf == NULL)
116 return(0); /* all diffuse */
117 sv = ndp->sd->rLambFront;
118 break;
119 case 0:
120 if (ndp->sd->rb == NULL)
121 return(0); /* all diffuse */
122 sv = ndp->sd->rLambBack;
123 break;
124 default:
125 if ((ndp->sd->tf == NULL) & (ndp->sd->tb == NULL))
126 return(0); /* all diffuse */
127 sv = ndp->sd->tLamb;
128 break;
129 }
130 if ((sv.cieY *= 1./PI) > FTINY) {
131 diffY = sv.cieY;
132 cvt_sdcolor(cdiff, &sv);
133 } else {
134 diffY = .0;
135 setcolor(cdiff, .0, .0, .0);
136 }
137 /* assign number of samples */
138 ec = SDsizeBSDF(&tomega, ndp->vray, vsrc, SDqueryMin, ndp->sd);
139 if (ec)
140 goto baderror;
141 /* check indirect over-counting */
142 if (ndp->thick != 0 && ndp->pr->crtype & (SPECULAR|AMBIENT)
143 && (vsrc[2] > 0) ^ (ndp->vray[2] > 0)) {
144 double dx = vsrc[0] + ndp->vray[0];
145 double dy = vsrc[1] + ndp->vray[1];
146 if (dx*dx + dy*dy <= omega+tomega)
147 return(0);
148 }
149 sf = specjitter * ndp->pr->rweight;
150 if (tomega <= .0)
151 nsamp = 1;
152 else if (25.*tomega <= omega)
153 nsamp = 100.*sf + .5;
154 else
155 nsamp = 4.*sf*omega/tomega + .5;
156 nsamp += !nsamp;
157 setcolor(cval, .0, .0, .0); /* sample our source area */
158 sf = sqrt(omega);
159 tsr = sqrt(tomega);
160 for (i = nsamp; i--; ) {
161 VCOPY(vsmp, vsrc); /* jitter query directions */
162 if (nsamp > 1) {
163 multisamp(sd, 2, (i + frandom())/(double)nsamp);
164 vsmp[0] += (sd[0] - .5)*sf;
165 vsmp[1] += (sd[1] - .5)*sf;
166 if (normalize(vsmp) == 0) {
167 --nsamp;
168 continue;
169 }
170 }
171 bsdf_jitter(vjit, ndp, tsr);
172 /* compute BSDF */
173 ec = SDevalBSDF(&sv, vjit, vsmp, ndp->sd);
174 if (ec)
175 goto baderror;
176 if (sv.cieY - diffY <= FTINY) {
177 addcolor(cval, cdiff);
178 continue; /* no specular part */
179 }
180 cvt_sdcolor(csmp, &sv);
181 addcolor(cval, csmp); /* else average it in */
182 ++ok;
183 }
184 if (!ok) {
185 setcolor(cval, .0, .0, .0);
186 return(0); /* no valid specular samples */
187 }
188 sf = 1./(double)nsamp;
189 scalecolor(cval, sf);
190 /* subtract diffuse contribution */
191 for (i = 3*(diffY > FTINY); i--; )
192 if ((colval(cval,i) -= colval(cdiff,i)) < .0)
193 colval(cval,i) = .0;
194 return(1);
195 baderror:
196 objerror(ndp->mp, USER, transSDError(ec));
197 return(0); /* gratis return */
198 }
199
200 /* Compute source contribution for BSDF (reflected & transmitted) */
201 static void
202 dir_bsdf(
203 COLOR cval, /* returned coefficient */
204 void *nnp, /* material data */
205 FVECT ldir, /* light source direction */
206 double omega /* light source size */
207 )
208 {
209 BSDFDAT *np = (BSDFDAT *)nnp;
210 double ldot;
211 double dtmp;
212 COLOR ctmp;
213
214 setcolor(cval, .0, .0, .0);
215
216 ldot = DOT(np->pnorm, ldir);
217 if ((-FTINY <= ldot) & (ldot <= FTINY))
218 return;
219
220 if (ldot > 0 && bright(np->rdiff) > FTINY) {
221 /*
222 * Compute added diffuse reflected component.
223 */
224 copycolor(ctmp, np->rdiff);
225 dtmp = ldot * omega * (1./PI);
226 scalecolor(ctmp, dtmp);
227 addcolor(cval, ctmp);
228 }
229 if (ldot < 0 && bright(np->tdiff) > FTINY) {
230 /*
231 * Compute added diffuse transmission.
232 */
233 copycolor(ctmp, np->tdiff);
234 dtmp = -ldot * omega * (1.0/PI);
235 scalecolor(ctmp, dtmp);
236 addcolor(cval, ctmp);
237 }
238 if (ambRayInPmap(np->pr))
239 return; /* specular already in photon map */
240 /*
241 * Compute scattering coefficient using BSDF.
242 */
243 if (!direct_bsdf_OK(ctmp, ldir, omega, np))
244 return;
245 if (ldot < 0) { /* pattern for specular transmission */
246 multcolor(ctmp, np->pr->pcol);
247 dtmp = -ldot * omega;
248 } else
249 dtmp = ldot * omega;
250 scalecolor(ctmp, dtmp);
251 addcolor(cval, ctmp);
252 }
253
254 /* Compute source contribution for BSDF (reflected only) */
255 static void
256 dir_brdf(
257 COLOR cval, /* returned coefficient */
258 void *nnp, /* material data */
259 FVECT ldir, /* light source direction */
260 double omega /* light source size */
261 )
262 {
263 BSDFDAT *np = (BSDFDAT *)nnp;
264 double ldot;
265 double dtmp;
266 COLOR ctmp, ctmp1, ctmp2;
267
268 setcolor(cval, .0, .0, .0);
269
270 ldot = DOT(np->pnorm, ldir);
271
272 if (ldot <= FTINY)
273 return;
274
275 if (bright(np->rdiff) > FTINY) {
276 /*
277 * Compute added diffuse reflected component.
278 */
279 copycolor(ctmp, np->rdiff);
280 dtmp = ldot * omega * (1./PI);
281 scalecolor(ctmp, dtmp);
282 addcolor(cval, ctmp);
283 }
284 if (ambRayInPmap(np->pr))
285 return; /* specular already in photon map */
286 /*
287 * Compute reflection coefficient using BSDF.
288 */
289 if (!direct_bsdf_OK(ctmp, ldir, omega, np))
290 return;
291 dtmp = ldot * omega;
292 scalecolor(ctmp, dtmp);
293 addcolor(cval, ctmp);
294 }
295
296 /* Compute source contribution for BSDF (transmitted only) */
297 static void
298 dir_btdf(
299 COLOR cval, /* returned coefficient */
300 void *nnp, /* material data */
301 FVECT ldir, /* light source direction */
302 double omega /* light source size */
303 )
304 {
305 BSDFDAT *np = (BSDFDAT *)nnp;
306 double ldot;
307 double dtmp;
308 COLOR ctmp;
309
310 setcolor(cval, .0, .0, .0);
311
312 ldot = DOT(np->pnorm, ldir);
313
314 if (ldot >= -FTINY)
315 return;
316
317 if (bright(np->tdiff) > FTINY) {
318 /*
319 * Compute added diffuse transmission.
320 */
321 copycolor(ctmp, np->tdiff);
322 dtmp = -ldot * omega * (1.0/PI);
323 scalecolor(ctmp, dtmp);
324 addcolor(cval, ctmp);
325 }
326 if (ambRayInPmap(np->pr))
327 return; /* specular already in photon map */
328 /*
329 * Compute scattering coefficient using BSDF.
330 */
331 if (!direct_bsdf_OK(ctmp, ldir, omega, np))
332 return;
333 /* full pattern on transmission */
334 multcolor(ctmp, np->pr->pcol);
335 dtmp = -ldot * omega;
336 scalecolor(ctmp, dtmp);
337 addcolor(cval, ctmp);
338 }
339
340 /* Sample separate BSDF component */
341 static int
342 sample_sdcomp(BSDFDAT *ndp, SDComponent *dcp, int usepat)
343 {
344 int nstarget = 1;
345 int nsent;
346 SDError ec;
347 SDValue bsv;
348 double xrand;
349 FVECT vsmp;
350 RAY sr;
351 /* multiple samples? */
352 if (specjitter > 1.5) {
353 nstarget = specjitter*ndp->pr->rweight + .5;
354 nstarget += !nstarget;
355 }
356 /* run through our samples */
357 for (nsent = 0; nsent < nstarget; nsent++) {
358 if (nstarget == 1) { /* stratify random variable */
359 xrand = urand(ilhash(dimlist,ndims)+samplendx);
360 if (specjitter < 1.)
361 xrand = .5 + specjitter*(xrand-.5);
362 } else {
363 xrand = (nsent + frandom())/(double)nstarget;
364 }
365 SDerrorDetail[0] = '\0'; /* sample direction & coef. */
366 bsdf_jitter(vsmp, ndp, ndp->sr_vpsa[0]);
367 ec = SDsampComponent(&bsv, vsmp, xrand, dcp);
368 if (ec)
369 objerror(ndp->mp, USER, transSDError(ec));
370 if (bsv.cieY <= FTINY) /* zero component? */
371 break;
372 /* map vector to world */
373 if (SDmapDir(sr.rdir, ndp->fromloc, vsmp) != SDEnone)
374 break;
375 /* spawn a specular ray */
376 if (nstarget > 1)
377 bsv.cieY /= (double)nstarget;
378 cvt_sdcolor(sr.rcoef, &bsv); /* use sample color */
379 if (usepat) /* apply pattern? */
380 multcolor(sr.rcoef, ndp->pr->pcol);
381 if (rayorigin(&sr, SPECULAR, ndp->pr, sr.rcoef) < 0) {
382 if (maxdepth > 0)
383 break;
384 continue; /* Russian roulette victim */
385 }
386 /* need to offset origin? */
387 if (ndp->thick != 0 && (ndp->pr->rod > 0) ^ (vsmp[2] > 0))
388 VSUM(sr.rorg, sr.rorg, ndp->pr->ron, -ndp->thick);
389 rayvalue(&sr); /* send & evaluate sample */
390 multcolor(sr.rcol, sr.rcoef);
391 addcolor(ndp->pr->rcol, sr.rcol);
392 }
393 return(nsent);
394 }
395
396 /* Sample non-diffuse components of BSDF */
397 static int
398 sample_sdf(BSDFDAT *ndp, int sflags)
399 {
400 int n, ntotal = 0;
401 SDSpectralDF *dfp;
402 COLORV *unsc;
403
404 if (sflags == SDsampSpT) {
405 unsc = ndp->tdiff;
406 if (ndp->pr->rod > 0)
407 dfp = (ndp->sd->tf != NULL) ? ndp->sd->tf : ndp->sd->tb;
408 else
409 dfp = (ndp->sd->tb != NULL) ? ndp->sd->tb : ndp->sd->tf;
410 } else /* sflags == SDsampSpR */ {
411 unsc = ndp->rdiff;
412 if (ndp->pr->rod > 0)
413 dfp = ndp->sd->rf;
414 else
415 dfp = ndp->sd->rb;
416 }
417 if (dfp == NULL) /* no specular component? */
418 return(0);
419 /* below sampling threshold? */
420 if (dfp->maxHemi <= specthresh+FTINY) {
421 if (dfp->maxHemi > FTINY) { /* XXX no color from BSDF */
422 FVECT vjit;
423 double d;
424 COLOR ctmp;
425 bsdf_jitter(vjit, ndp, ndp->sr_vpsa[1]);
426 d = SDdirectHemi(vjit, sflags, ndp->sd);
427 if (sflags == SDsampSpT) {
428 copycolor(ctmp, ndp->pr->pcol);
429 scalecolor(ctmp, d);
430 } else /* no pattern on reflection */
431 setcolor(ctmp, d, d, d);
432 addcolor(unsc, ctmp);
433 }
434 return(0);
435 }
436 /* else need to sample */
437 dimlist[ndims++] = (int)(size_t)ndp->mp;
438 ndims++;
439 for (n = dfp->ncomp; n--; ) { /* loop over components */
440 dimlist[ndims-1] = n + 9438;
441 ntotal += sample_sdcomp(ndp, &dfp->comp[n], sflags==SDsampSpT);
442 }
443 ndims -= 2;
444 return(ntotal);
445 }
446
447 /* Color a ray that hit a BSDF material */
448 int
449 m_bsdf(OBJREC *m, RAY *r)
450 {
451 int hitfront;
452 COLOR ctmp;
453 SDError ec;
454 FVECT upvec, vtmp;
455 MFUNC *mf;
456 BSDFDAT nd;
457 /* check arguments */
458 if ((m->oargs.nsargs < 6) | (m->oargs.nfargs > 9) |
459 (m->oargs.nfargs % 3))
460 objerror(m, USER, "bad # arguments");
461 /* record surface struck */
462 hitfront = (r->rod > 0);
463 /* load cal file */
464 mf = getfunc(m, 5, 0x1d, 1);
465 setfunc(m, r);
466 /* get thickness */
467 nd.thick = evalue(mf->ep[0]);
468 if ((-FTINY <= nd.thick) & (nd.thick <= FTINY))
469 nd.thick = .0;
470 /* check shadow */
471 if (r->crtype & SHADOW) {
472 if (nd.thick != 0)
473 raytrans(r); /* pass-through */
474 return(1); /* or shadow */
475 }
476 /* check backface visibility */
477 if (!hitfront & !backvis) {
478 raytrans(r);
479 return(1);
480 }
481 /* check other rays to pass */
482 if (nd.thick != 0 && (!(r->crtype & (SPECULAR|AMBIENT)) ||
483 (nd.thick > 0) ^ hitfront)) {
484 raytrans(r); /* hide our proxy */
485 return(1);
486 }
487 nd.mp = m;
488 nd.pr = r;
489 /* get BSDF data */
490 nd.sd = loadBSDF(m->oargs.sarg[1]);
491 /* diffuse reflectance */
492 if (hitfront) {
493 cvt_sdcolor(nd.rdiff, &nd.sd->rLambFront);
494 if (m->oargs.nfargs >= 3) {
495 setcolor(ctmp, m->oargs.farg[0],
496 m->oargs.farg[1],
497 m->oargs.farg[2]);
498 addcolor(nd.rdiff, ctmp);
499 }
500 } else {
501 cvt_sdcolor(nd.rdiff, &nd.sd->rLambBack);
502 if (m->oargs.nfargs >= 6) {
503 setcolor(ctmp, m->oargs.farg[3],
504 m->oargs.farg[4],
505 m->oargs.farg[5]);
506 addcolor(nd.rdiff, ctmp);
507 }
508 }
509 /* diffuse transmittance */
510 cvt_sdcolor(nd.tdiff, &nd.sd->tLamb);
511 if (m->oargs.nfargs >= 9) {
512 setcolor(ctmp, m->oargs.farg[6],
513 m->oargs.farg[7],
514 m->oargs.farg[8]);
515 addcolor(nd.tdiff, ctmp);
516 }
517 /* get modifiers */
518 raytexture(r, m->omod);
519 /* modify diffuse values */
520 multcolor(nd.rdiff, r->pcol);
521 multcolor(nd.tdiff, r->pcol);
522 /* get up vector */
523 upvec[0] = evalue(mf->ep[1]);
524 upvec[1] = evalue(mf->ep[2]);
525 upvec[2] = evalue(mf->ep[3]);
526 /* return to world coords */
527 if (mf->fxp != &unitxf) {
528 multv3(upvec, upvec, mf->fxp->xfm);
529 nd.thick *= mf->fxp->sca;
530 }
531 if (r->rox != NULL) {
532 multv3(upvec, upvec, r->rox->f.xfm);
533 nd.thick *= r->rox->f.sca;
534 }
535 raynormal(nd.pnorm, r);
536 /* compute local BSDF xform */
537 ec = SDcompXform(nd.toloc, nd.pnorm, upvec);
538 if (!ec) {
539 nd.vray[0] = -r->rdir[0];
540 nd.vray[1] = -r->rdir[1];
541 nd.vray[2] = -r->rdir[2];
542 ec = SDmapDir(nd.vray, nd.toloc, nd.vray);
543 }
544 if (!ec)
545 ec = SDinvXform(nd.fromloc, nd.toloc);
546 if (ec) {
547 objerror(m, WARNING, "Illegal orientation vector");
548 return(1);
549 }
550 /* determine BSDF resolution */
551 ec = SDsizeBSDF(nd.sr_vpsa, nd.vray, NULL, SDqueryMin+SDqueryMax, nd.sd);
552 if (ec)
553 objerror(m, USER, transSDError(ec));
554
555 nd.sr_vpsa[0] = sqrt(nd.sr_vpsa[0]);
556 nd.sr_vpsa[1] = sqrt(nd.sr_vpsa[1]);
557 if (!hitfront) { /* perturb normal towards hit */
558 nd.pnorm[0] = -nd.pnorm[0];
559 nd.pnorm[1] = -nd.pnorm[1];
560 nd.pnorm[2] = -nd.pnorm[2];
561 }
562 /* sample reflection */
563 sample_sdf(&nd, SDsampSpR);
564 /* sample transmission */
565 sample_sdf(&nd, SDsampSpT);
566 /* compute indirect diffuse */
567 if (bright(nd.rdiff) > FTINY) { /* ambient from reflection */
568 if (!hitfront)
569 flipsurface(r);
570 copycolor(ctmp, nd.rdiff);
571 multambient(ctmp, r, nd.pnorm);
572 addcolor(r->rcol, ctmp);
573 if (!hitfront)
574 flipsurface(r);
575 }
576 if (bright(nd.tdiff) > FTINY) { /* ambient from other side */
577 FVECT bnorm;
578 if (hitfront)
579 flipsurface(r);
580 bnorm[0] = -nd.pnorm[0];
581 bnorm[1] = -nd.pnorm[1];
582 bnorm[2] = -nd.pnorm[2];
583 copycolor(ctmp, nd.tdiff);
584 if (nd.thick != 0) { /* proxy with offset? */
585 VCOPY(vtmp, r->rop);
586 VSUM(r->rop, vtmp, r->ron, nd.thick);
587 multambient(ctmp, r, bnorm);
588 VCOPY(r->rop, vtmp);
589 } else
590 multambient(ctmp, r, bnorm);
591 addcolor(r->rcol, ctmp);
592 if (hitfront)
593 flipsurface(r);
594 }
595 /* add direct component */
596 if ((bright(nd.tdiff) <= FTINY) & (nd.sd->tf == NULL) &
597 (nd.sd->tb == NULL)) {
598 direct(r, dir_brdf, &nd); /* reflection only */
599 } else if (nd.thick == 0) {
600 direct(r, dir_bsdf, &nd); /* thin surface scattering */
601 } else {
602 direct(r, dir_brdf, &nd); /* reflection first */
603 VCOPY(vtmp, r->rop); /* offset for transmitted */
604 VSUM(r->rop, vtmp, r->ron, -nd.thick);
605 direct(r, dir_btdf, &nd); /* separate transmission */
606 VCOPY(r->rop, vtmp);
607 }
608 /* clean up */
609 SDfreeCache(nd.sd);
610 return(1);
611 }