ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/rt/m_bsdf.c
Revision: 2.59
Committed: Wed Jun 3 02:27:32 2020 UTC (4 years ago) by greg
Content type: text/plain
Branch: MAIN
Changes since 2.58: +5 -4 lines
Log Message:
Removed small source sampling bias in direct_specular_OK()

File Contents

# User Rev Content
1 greg 2.1 #ifndef lint
2 greg 2.59 static const char RCSid[] = "$Id: m_bsdf.c,v 2.58 2020/06/03 02:08:32 greg Exp $";
3 greg 2.1 #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 greg 2.50 #include "otypes.h"
12 greg 2.1 #include "ambient.h"
13     #include "source.h"
14     #include "func.h"
15     #include "bsdf.h"
16     #include "random.h"
17 greg 2.30 #include "pmapmat.h"
18 greg 2.1
19     /*
20 greg 2.50 * Arguments to this material include optional diffuse colors.
21 greg 2.1 * String arguments include the BSDF and function files.
22 greg 2.50 * For the MAT_BSDF type, a non-zero thickness causes the useful behavior
23 greg 2.5 * of translating transmitted rays this distance beneath the surface
24     * (opposite the surface normal) to bypass any intervening geometry.
25     * Translation only affects scattered, non-source-directed samples.
26     * A non-zero thickness has the further side-effect that an unscattered
27 greg 2.35 * (view) ray will pass right through our material, making the BSDF
28     * surface invisible and showing the proxied geometry instead. Thickness
29     * has the further effect of turning off reflection on the reverse side so
30     * rays heading in the opposite direction pass unimpeded through the BSDF
31 greg 2.5 * 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 greg 2.52 * For the MAT_ABSDF type, we check for a strong "through" component.
40 greg 2.50 * Such a component will cause direct rays to pass through unscattered.
41 greg 2.40 * A separate test prevents over-counting by dropping samples that are
42     * too close to this "through" direction. BSDFs with such a through direction
43     * will also have a view component, meaning they are somewhat see-through.
44 greg 2.52 * A MAT_BSDF type with zero thickness behaves the same as a MAT_ABSDF
45 greg 2.50 * type with no strong through component.
46 greg 2.1 * The "up" vector for the BSDF is given by three variables, defined
47     * (along with the thickness) by the named function file, or '.' if none.
48     * Together with the surface normal, this defines the local coordinate
49     * system for the BSDF.
50     * We do not reorient the surface, so if the BSDF has no back-side
51 greg 2.5 * reflectance and none is given in the real arguments, a BSDF surface
52     * with zero thickness will appear black when viewed from behind
53 greg 2.35 * unless backface visibility is on, when it becomes invisible.
54 greg 2.5 * The diffuse arguments are added to components in the BSDF file,
55 greg 2.1 * not multiplied. However, patterns affect this material as a multiplier
56     * on everything except non-diffuse reflection.
57     *
58 greg 2.52 * Arguments for MAT_ABSDF are:
59 greg 2.50 * 5+ BSDFfile ux uy uz funcfile transform
60     * 0
61     * 0|3|6|9 rdf gdf bdf
62     * rdb gdb bdb
63     * rdt gdt bdt
64     *
65 greg 2.1 * Arguments for MAT_BSDF are:
66     * 6+ thick BSDFfile ux uy uz funcfile transform
67     * 0
68 greg 2.8 * 0|3|6|9 rdf gdf bdf
69 greg 2.1 * rdb gdb bdb
70     * rdt gdt bdt
71     */
72    
73 greg 2.4 /*
74     * Note that our reverse ray-tracing process means that the positions
75     * of incoming and outgoing vectors may be reversed in our calls
76 greg 2.35 * to the BSDF library. This is usually fine, since the bidirectional nature
77 greg 2.4 * of the BSDF (that's what the 'B' stands for) means it all works out.
78     */
79    
80 greg 2.1 typedef struct {
81     OBJREC *mp; /* material pointer */
82     RAY *pr; /* intersected ray */
83     FVECT pnorm; /* perturbed surface normal */
84 greg 2.4 FVECT vray; /* local outgoing (return) vector */
85 greg 2.9 double sr_vpsa[2]; /* sqrt of BSDF projected solid angle extrema */
86 greg 2.1 RREAL toloc[3][3]; /* world to local BSDF coords */
87     RREAL fromloc[3][3]; /* local BSDF coords to world */
88     double thick; /* surface thickness */
89 greg 2.52 COLOR cthru; /* "through" component for MAT_ABSDF */
90 greg 2.1 SDData *sd; /* loaded BSDF data */
91 greg 2.31 COLOR rdiff; /* diffuse reflection */
92 greg 2.39 COLOR runsamp; /* BSDF hemispherical reflection */
93 greg 2.31 COLOR tdiff; /* diffuse transmission */
94 greg 2.39 COLOR tunsamp; /* BSDF hemispherical transmission */
95 greg 2.1 } BSDFDAT; /* BSDF material data */
96    
97     #define cvt_sdcolor(cv, svp) ccy2rgb(&(svp)->spec, (svp)->cieY, cv)
98    
99 greg 2.58 typedef struct {
100     double vy; /* brightness (for sorting) */
101     FVECT tdir; /* through sample direction (normalized) */
102     COLOR vcol; /* BTDF color */
103     } PEAKSAMP; /* BTDF peak sample */
104    
105     /* Comparison function to put near-peak values in descending order */
106     static int
107     cmp_psamp(const void *p1, const void *p2)
108     {
109     double diff = (*(const PEAKSAMP *)p1).vy - (*(const PEAKSAMP *)p2).vy;
110     if (diff > 0) return(-1);
111     if (diff < 0) return(1);
112     return(0);
113     }
114    
115 greg 2.52 /* Compute "through" component color for MAT_ABSDF */
116 greg 2.34 static void
117     compute_through(BSDFDAT *ndp)
118     {
119     #define NDIR2CHECK 13
120     static const float dir2check[NDIR2CHECK][2] = {
121     {0, 0},
122     {-0.8, 0},
123     {0, 0.8},
124     {0, -0.8},
125     {0.8, 0},
126     {-0.8, 0.8},
127     {-0.8, -0.8},
128     {0.8, 0.8},
129     {0.8, -0.8},
130     {-1.6, 0},
131     {0, 1.6},
132     {0, -1.6},
133     {1.6, 0},
134     };
135 greg 2.58 const double peak_over = 1.5;
136     PEAKSAMP psamp[NDIR2CHECK];
137 greg 2.34 SDSpectralDF *dfp;
138     FVECT pdir;
139     double tomega, srchrad;
140 greg 2.58 double tomsum;
141     COLOR vpeak;
142     double vypeak, vysum;
143     int i, ns, ntot;
144 greg 2.34 SDError ec;
145    
146     if (ndp->pr->rod > 0)
147     dfp = (ndp->sd->tf != NULL) ? ndp->sd->tf : ndp->sd->tb;
148     else
149     dfp = (ndp->sd->tb != NULL) ? ndp->sd->tb : ndp->sd->tf;
150    
151     if (dfp == NULL)
152     return; /* no specular transmission */
153     if (bright(ndp->pr->pcol) <= FTINY)
154     return; /* pattern is black, here */
155 greg 2.58 srchrad = sqrt(dfp->minProjSA); /* else evaluate peak */
156     vysum = 0;
157     for (i = 0; i < NDIR2CHECK; i++) {
158     SDValue sv;
159     psamp[i].tdir[0] = -ndp->vray[0] + dir2check[i][0]*srchrad;
160     psamp[i].tdir[1] = -ndp->vray[1] + dir2check[i][1]*srchrad;
161     psamp[i].tdir[2] = -ndp->vray[2];
162     normalize(psamp[i].tdir);
163     ec = SDevalBSDF(&sv, psamp[i].tdir, ndp->vray, ndp->sd);
164     if (ec)
165     goto baderror;
166     cvt_sdcolor(psamp[i].vcol, &sv);
167     vysum += psamp[i].vy = sv.cieY;
168     }
169     if (vysum <= FTINY) /* zero neighborhood? */
170     return;
171     qsort(psamp, NDIR2CHECK, sizeof(PEAKSAMP), cmp_psamp);
172 greg 2.40 setcolor(vpeak, 0, 0, 0);
173 greg 2.58 vypeak = tomsum = 0; /* combine top unique values */
174     ns = 0; ntot = NDIR2CHECK;
175 greg 2.34 for (i = 0; i < NDIR2CHECK; i++) {
176 greg 2.58 if (i) {
177     if (psamp[i].vy == psamp[i-1].vy) {
178     vysum -= psamp[i].vy;
179     --ntot;
180     continue; /* assume duplicate sample */
181     }
182     if (vypeak > 8.*psamp[i].vy*ns)
183     continue; /* peak cut-off */
184     }
185     ec = SDsizeBSDF(&tomega, psamp[i].tdir, ndp->vray,
186     SDqueryMin, ndp->sd);
187 greg 2.34 if (ec)
188     goto baderror;
189 greg 2.58 if (tomega > 1.5*dfp->minProjSA) {
190     if (!i) return; /* not really a peak? */
191     continue;
192 greg 2.34 }
193 greg 2.58 scalecolor(psamp[i].vcol, tomega);
194     addcolor(vpeak, psamp[i].vcol);
195     tomsum += tomega;
196     vypeak += psamp[i].vy;
197     ++ns;
198     }
199     if (vypeak*(ntot-ns) < peak_over*(vysum-vypeak)*ns)
200     return; /* peak not peaky enough */
201     if ((vypeak/ns - ndp->sd->tLamb.cieY*(1./PI))*tomsum <= .001)
202 greg 2.40 return; /* < 0.1% transmission */
203 greg 2.58 copycolor(ndp->cthru, vpeak); /* already scaled by omega */
204 greg 2.34 multcolor(ndp->cthru, ndp->pr->pcol); /* modify by pattern */
205     return;
206     baderror:
207     objerror(ndp->mp, USER, transSDError(ec));
208     #undef NDIR2CHECK
209     }
210    
211 greg 2.4 /* Jitter ray sample according to projected solid angle and specjitter */
212     static void
213 greg 2.15 bsdf_jitter(FVECT vres, BSDFDAT *ndp, double sr_psa)
214 greg 2.4 {
215     VCOPY(vres, ndp->vray);
216     if (specjitter < 1.)
217     sr_psa *= specjitter;
218     if (sr_psa <= FTINY)
219     return;
220     vres[0] += sr_psa*(.5 - frandom());
221     vres[1] += sr_psa*(.5 - frandom());
222     normalize(vres);
223     }
224    
225 greg 2.33 /* Get BSDF specular for direct component, returning true if OK to proceed */
226 greg 2.7 static int
227 greg 2.33 direct_specular_OK(COLOR cval, FVECT ldir, double omega, BSDFDAT *ndp)
228 greg 2.7 {
229 greg 2.43 int nsamp;
230     double wtot = 0;
231 greg 2.13 FVECT vsrc, vsmp, vjit;
232 greg 2.36 double tomega, tomega2;
233 greg 2.15 double sf, tsr, sd[2];
234 greg 2.32 COLOR csmp, cdiff;
235     double diffY;
236 greg 2.7 SDValue sv;
237     SDError ec;
238 greg 2.13 int i;
239 greg 2.37 /* in case we fail */
240 greg 2.40 setcolor(cval, 0, 0, 0);
241 greg 2.7 /* transform source direction */
242     if (SDmapDir(vsrc, ndp->toloc, ldir) != SDEnone)
243     return(0);
244 greg 2.32 /* will discount diffuse portion */
245     switch ((vsrc[2] > 0)<<1 | (ndp->vray[2] > 0)) {
246     case 3:
247     if (ndp->sd->rf == NULL)
248     return(0); /* all diffuse */
249     sv = ndp->sd->rLambFront;
250     break;
251     case 0:
252     if (ndp->sd->rb == NULL)
253     return(0); /* all diffuse */
254     sv = ndp->sd->rLambBack;
255     break;
256     default:
257     if ((ndp->sd->tf == NULL) & (ndp->sd->tb == NULL))
258     return(0); /* all diffuse */
259     sv = ndp->sd->tLamb;
260     break;
261     }
262 greg 2.33 if (sv.cieY > FTINY) {
263     diffY = sv.cieY *= 1./PI;
264 greg 2.32 cvt_sdcolor(cdiff, &sv);
265     } else {
266 greg 2.40 diffY = 0;
267     setcolor(cdiff, 0, 0, 0);
268 greg 2.32 }
269 greg 2.53 /* need projected solid angle */
270 greg 2.37 omega *= fabs(vsrc[2]);
271 greg 2.13 /* check indirect over-counting */
272 greg 2.40 if ((vsrc[2] > 0) ^ (ndp->vray[2] > 0) && bright(ndp->cthru) > FTINY) {
273 greg 2.53 double dx = vsrc[0] + ndp->vray[0];
274     double dy = vsrc[1] + ndp->vray[1];
275     SDSpectralDF *dfp = (ndp->pr->rod > 0) ?
276     ((ndp->sd->tf != NULL) ? ndp->sd->tf : ndp->sd->tb) :
277     ((ndp->sd->tb != NULL) ? ndp->sd->tb : ndp->sd->tf) ;
278    
279     if (dx*dx + dy*dy <= (2.5*4./PI)*(omega + dfp->minProjSA +
280     2.*sqrt(omega*dfp->minProjSA)))
281 greg 2.7 return(0);
282     }
283 greg 2.53 ec = SDsizeBSDF(&tomega, ndp->vray, vsrc, SDqueryMin, ndp->sd);
284     if (ec)
285     goto baderror;
286 greg 2.37 /* assign number of samples */
287 greg 2.15 sf = specjitter * ndp->pr->rweight;
288 greg 2.40 if (tomega <= 0)
289 greg 2.24 nsamp = 1;
290     else if (25.*tomega <= omega)
291 greg 2.15 nsamp = 100.*sf + .5;
292     else
293     nsamp = 4.*sf*omega/tomega + .5;
294     nsamp += !nsamp;
295 greg 2.37 sf = sqrt(omega); /* sample our source area */
296 greg 2.15 tsr = sqrt(tomega);
297 greg 2.13 for (i = nsamp; i--; ) {
298     VCOPY(vsmp, vsrc); /* jitter query directions */
299     if (nsamp > 1) {
300     multisamp(sd, 2, (i + frandom())/(double)nsamp);
301     vsmp[0] += (sd[0] - .5)*sf;
302     vsmp[1] += (sd[1] - .5)*sf;
303 greg 2.36 normalize(vsmp);
304 greg 2.13 }
305 greg 2.15 bsdf_jitter(vjit, ndp, tsr);
306 greg 2.37 /* compute BSDF */
307     ec = SDevalBSDF(&sv, vjit, vsmp, ndp->sd);
308     if (ec)
309     goto baderror;
310     if (sv.cieY - diffY <= FTINY)
311     continue; /* no specular part */
312 greg 2.36 /* check for variable resolution */
313     ec = SDsizeBSDF(&tomega2, vjit, vsmp, SDqueryMin, ndp->sd);
314     if (ec)
315     goto baderror;
316     if (tomega2 < .12*tomega)
317     continue; /* not safe to include */
318 greg 2.13 cvt_sdcolor(csmp, &sv);
319 greg 2.59 #if 0
320     if (sf < 2.5*tsr) { /* weight by BSDF for small sources */
321 greg 2.44 scalecolor(csmp, sv.cieY);
322     wtot += sv.cieY;
323     } else
324 greg 2.59 #endif
325     wtot += 1.;
326 greg 2.43 addcolor(cval, csmp);
327 greg 2.13 }
328 greg 2.43 if (wtot <= FTINY) /* no valid specular samples? */
329 greg 2.37 return(0);
330    
331 greg 2.43 sf = 1./wtot; /* weighted average BSDF */
332 greg 2.13 scalecolor(cval, sf);
333 greg 2.32 /* subtract diffuse contribution */
334     for (i = 3*(diffY > FTINY); i--; )
335 greg 2.40 if ((colval(cval,i) -= colval(cdiff,i)) < 0)
336     colval(cval,i) = 0;
337 greg 2.32 return(1);
338 greg 2.13 baderror:
339     objerror(ndp->mp, USER, transSDError(ec));
340 greg 2.17 return(0); /* gratis return */
341 greg 2.7 }
342    
343 greg 2.5 /* Compute source contribution for BSDF (reflected & transmitted) */
344 greg 2.1 static void
345 greg 2.5 dir_bsdf(
346 greg 2.1 COLOR cval, /* returned coefficient */
347     void *nnp, /* material data */
348     FVECT ldir, /* light source direction */
349     double omega /* light source size */
350     )
351     {
352 greg 2.3 BSDFDAT *np = (BSDFDAT *)nnp;
353 greg 2.1 double ldot;
354     double dtmp;
355     COLOR ctmp;
356    
357 greg 2.40 setcolor(cval, 0, 0, 0);
358 greg 2.1
359     ldot = DOT(np->pnorm, ldir);
360     if ((-FTINY <= ldot) & (ldot <= FTINY))
361     return;
362    
363 greg 2.9 if (ldot > 0 && bright(np->rdiff) > FTINY) {
364 greg 2.1 /*
365 greg 2.39 * Compute diffuse reflected component
366 greg 2.1 */
367     copycolor(ctmp, np->rdiff);
368     dtmp = ldot * omega * (1./PI);
369     scalecolor(ctmp, dtmp);
370     addcolor(cval, ctmp);
371     }
372 greg 2.9 if (ldot < 0 && bright(np->tdiff) > FTINY) {
373 greg 2.1 /*
374 greg 2.39 * Compute diffuse transmission
375 greg 2.1 */
376     copycolor(ctmp, np->tdiff);
377     dtmp = -ldot * omega * (1.0/PI);
378     scalecolor(ctmp, dtmp);
379     addcolor(cval, ctmp);
380     }
381 greg 2.30 if (ambRayInPmap(np->pr))
382     return; /* specular already in photon map */
383 greg 2.1 /*
384 greg 2.39 * Compute specular scattering coefficient using BSDF
385 greg 2.1 */
386 greg 2.33 if (!direct_specular_OK(ctmp, ldir, omega, np))
387 greg 2.1 return;
388 greg 2.31 if (ldot < 0) { /* pattern for specular transmission */
389 greg 2.1 multcolor(ctmp, np->pr->pcol);
390     dtmp = -ldot * omega;
391 greg 2.31 } else
392     dtmp = ldot * omega;
393 greg 2.1 scalecolor(ctmp, dtmp);
394     addcolor(cval, ctmp);
395     }
396    
397 greg 2.5 /* Compute source contribution for BSDF (reflected only) */
398     static void
399     dir_brdf(
400     COLOR cval, /* returned coefficient */
401     void *nnp, /* material data */
402     FVECT ldir, /* light source direction */
403     double omega /* light source size */
404     )
405     {
406     BSDFDAT *np = (BSDFDAT *)nnp;
407     double ldot;
408     double dtmp;
409     COLOR ctmp, ctmp1, ctmp2;
410    
411 greg 2.40 setcolor(cval, 0, 0, 0);
412 greg 2.5
413     ldot = DOT(np->pnorm, ldir);
414    
415     if (ldot <= FTINY)
416     return;
417    
418     if (bright(np->rdiff) > FTINY) {
419     /*
420 greg 2.39 * Compute diffuse reflected component
421 greg 2.5 */
422     copycolor(ctmp, np->rdiff);
423     dtmp = ldot * omega * (1./PI);
424     scalecolor(ctmp, dtmp);
425     addcolor(cval, ctmp);
426     }
427 greg 2.30 if (ambRayInPmap(np->pr))
428     return; /* specular already in photon map */
429 greg 2.5 /*
430 greg 2.39 * Compute specular reflection coefficient using BSDF
431 greg 2.5 */
432 greg 2.33 if (!direct_specular_OK(ctmp, ldir, omega, np))
433 greg 2.5 return;
434     dtmp = ldot * omega;
435     scalecolor(ctmp, dtmp);
436     addcolor(cval, ctmp);
437     }
438    
439     /* Compute source contribution for BSDF (transmitted only) */
440     static void
441     dir_btdf(
442     COLOR cval, /* returned coefficient */
443     void *nnp, /* material data */
444     FVECT ldir, /* light source direction */
445     double omega /* light source size */
446     )
447     {
448     BSDFDAT *np = (BSDFDAT *)nnp;
449     double ldot;
450     double dtmp;
451     COLOR ctmp;
452    
453 greg 2.40 setcolor(cval, 0, 0, 0);
454 greg 2.5
455     ldot = DOT(np->pnorm, ldir);
456    
457     if (ldot >= -FTINY)
458     return;
459    
460     if (bright(np->tdiff) > FTINY) {
461     /*
462 greg 2.39 * Compute diffuse transmission
463 greg 2.5 */
464     copycolor(ctmp, np->tdiff);
465     dtmp = -ldot * omega * (1.0/PI);
466     scalecolor(ctmp, dtmp);
467     addcolor(cval, ctmp);
468     }
469 greg 2.30 if (ambRayInPmap(np->pr))
470     return; /* specular already in photon map */
471 greg 2.5 /*
472 greg 2.39 * Compute specular scattering coefficient using BSDF
473 greg 2.5 */
474 greg 2.33 if (!direct_specular_OK(ctmp, ldir, omega, np))
475 greg 2.5 return;
476     /* full pattern on transmission */
477     multcolor(ctmp, np->pr->pcol);
478     dtmp = -ldot * omega;
479     scalecolor(ctmp, dtmp);
480     addcolor(cval, ctmp);
481     }
482    
483 greg 2.1 /* Sample separate BSDF component */
484     static int
485 greg 2.40 sample_sdcomp(BSDFDAT *ndp, SDComponent *dcp, int xmit)
486 greg 2.1 {
487 greg 2.47 const int hasthru = (xmit &&
488     !(ndp->pr->crtype & (SPECULAR|AMBIENT))
489     && bright(ndp->cthru) > FTINY);
490 greg 2.41 int nstarget = 1;
491     int nsent = 0;
492     int n;
493     SDError ec;
494     SDValue bsv;
495     double xrand;
496     FVECT vsmp, vinc;
497     RAY sr;
498 greg 2.1 /* multiple samples? */
499     if (specjitter > 1.5) {
500     nstarget = specjitter*ndp->pr->rweight + .5;
501 greg 2.14 nstarget += !nstarget;
502 greg 2.1 }
503 greg 2.11 /* run through our samples */
504 greg 2.40 for (n = 0; n < nstarget; n++) {
505 greg 2.15 if (nstarget == 1) { /* stratify random variable */
506 greg 2.11 xrand = urand(ilhash(dimlist,ndims)+samplendx);
507 greg 2.15 if (specjitter < 1.)
508     xrand = .5 + specjitter*(xrand-.5);
509     } else {
510 greg 2.40 xrand = (n + frandom())/(double)nstarget;
511 greg 2.15 }
512 greg 2.11 SDerrorDetail[0] = '\0'; /* sample direction & coef. */
513 greg 2.15 bsdf_jitter(vsmp, ndp, ndp->sr_vpsa[0]);
514 greg 2.40 VCOPY(vinc, vsmp); /* to compare after */
515 greg 2.11 ec = SDsampComponent(&bsv, vsmp, xrand, dcp);
516 greg 2.1 if (ec)
517 greg 2.2 objerror(ndp->mp, USER, transSDError(ec));
518 greg 2.11 if (bsv.cieY <= FTINY) /* zero component? */
519 greg 2.1 break;
520 greg 2.40 if (hasthru) { /* check for view ray */
521     double dx = vinc[0] + vsmp[0];
522     double dy = vinc[1] + vsmp[1];
523     if (dx*dx + dy*dy <= ndp->sr_vpsa[0]*ndp->sr_vpsa[0])
524     continue; /* exclude view sample */
525     }
526     /* map non-view sample->world */
527 greg 2.4 if (SDmapDir(sr.rdir, ndp->fromloc, vsmp) != SDEnone)
528 greg 2.1 break;
529     /* spawn a specular ray */
530     if (nstarget > 1)
531     bsv.cieY /= (double)nstarget;
532 greg 2.11 cvt_sdcolor(sr.rcoef, &bsv); /* use sample color */
533 greg 2.40 if (xmit) /* apply pattern on transmit */
534 greg 2.1 multcolor(sr.rcoef, ndp->pr->pcol);
535     if (rayorigin(&sr, SPECULAR, ndp->pr, sr.rcoef) < 0) {
536 greg 2.48 if (!n & (nstarget > 1)) {
537 greg 2.49 n = nstarget; /* avoid infinitue loop */
538 greg 2.48 nstarget = nstarget*sr.rweight/minweight;
539 greg 2.49 if (n == nstarget) break;
540 greg 2.48 n = -1; /* moved target */
541     }
542     continue; /* try again */
543 greg 2.1 }
544 greg 2.40 if (xmit && ndp->thick != 0) /* need to offset origin? */
545 greg 2.5 VSUM(sr.rorg, sr.rorg, ndp->pr->ron, -ndp->thick);
546 greg 2.1 rayvalue(&sr); /* send & evaluate sample */
547     multcolor(sr.rcol, sr.rcoef);
548     addcolor(ndp->pr->rcol, sr.rcol);
549 greg 2.40 ++nsent;
550 greg 2.1 }
551     return(nsent);
552     }
553    
554     /* Sample non-diffuse components of BSDF */
555     static int
556     sample_sdf(BSDFDAT *ndp, int sflags)
557     {
558 greg 2.46 int hasthru = (sflags == SDsampSpT &&
559 greg 2.47 !(ndp->pr->crtype & (SPECULAR|AMBIENT))
560     && bright(ndp->cthru) > FTINY);
561 greg 2.1 int n, ntotal = 0;
562 greg 2.40 double b = 0;
563 greg 2.1 SDSpectralDF *dfp;
564     COLORV *unsc;
565    
566     if (sflags == SDsampSpT) {
567 greg 2.39 unsc = ndp->tunsamp;
568 greg 2.22 if (ndp->pr->rod > 0)
569     dfp = (ndp->sd->tf != NULL) ? ndp->sd->tf : ndp->sd->tb;
570     else
571     dfp = (ndp->sd->tb != NULL) ? ndp->sd->tb : ndp->sd->tf;
572 greg 2.1 } else /* sflags == SDsampSpR */ {
573 greg 2.39 unsc = ndp->runsamp;
574 greg 2.31 if (ndp->pr->rod > 0)
575 greg 2.1 dfp = ndp->sd->rf;
576 greg 2.31 else
577 greg 2.1 dfp = ndp->sd->rb;
578     }
579 greg 2.40 setcolor(unsc, 0, 0, 0);
580 greg 2.1 if (dfp == NULL) /* no specular component? */
581     return(0);
582 greg 2.40
583     if (hasthru) { /* separate view sample? */
584     RAY tr;
585     if (rayorigin(&tr, TRANS, ndp->pr, ndp->cthru) == 0) {
586     VCOPY(tr.rdir, ndp->pr->rdir);
587     rayvalue(&tr);
588     multcolor(tr.rcol, tr.rcoef);
589     addcolor(ndp->pr->rcol, tr.rcol);
590 greg 2.56 ndp->pr->rxt = ndp->pr->rot + raydistance(&tr);
591 greg 2.40 ++ntotal;
592     b = bright(ndp->cthru);
593     } else
594     hasthru = 0;
595     }
596 greg 2.43 if (dfp->maxHemi - b <= FTINY) { /* have specular to sample? */
597 greg 2.40 b = 0;
598     } else {
599     FVECT vjit;
600     bsdf_jitter(vjit, ndp, ndp->sr_vpsa[1]);
601     b = SDdirectHemi(vjit, sflags, ndp->sd) - b;
602     if (b < 0) b = 0;
603     }
604     if (b <= specthresh+FTINY) { /* below sampling threshold? */
605     if (b > FTINY) { /* XXX no color from BSDF */
606 greg 2.1 if (sflags == SDsampSpT) {
607 greg 2.39 copycolor(unsc, ndp->pr->pcol);
608 greg 2.40 scalecolor(unsc, b);
609 greg 2.1 } else /* no pattern on reflection */
610 greg 2.40 setcolor(unsc, b, b, b);
611 greg 2.1 }
612 greg 2.40 return(ntotal);
613 greg 2.1 }
614 greg 2.41 dimlist[ndims] = (int)(size_t)ndp->mp; /* else sample specular */
615     ndims += 2;
616 greg 2.1 for (n = dfp->ncomp; n--; ) { /* loop over components */
617     dimlist[ndims-1] = n + 9438;
618     ntotal += sample_sdcomp(ndp, &dfp->comp[n], sflags==SDsampSpT);
619     }
620     ndims -= 2;
621     return(ntotal);
622     }
623    
624     /* Color a ray that hit a BSDF material */
625     int
626     m_bsdf(OBJREC *m, RAY *r)
627     {
628 greg 2.50 int hasthick = (m->otype == MAT_BSDF);
629 greg 2.6 int hitfront;
630 greg 2.1 COLOR ctmp;
631     SDError ec;
632 greg 2.5 FVECT upvec, vtmp;
633 greg 2.1 MFUNC *mf;
634     BSDFDAT nd;
635     /* check arguments */
636 greg 2.50 if ((m->oargs.nsargs < hasthick+5) | (m->oargs.nfargs > 9) |
637 greg 2.1 (m->oargs.nfargs % 3))
638     objerror(m, USER, "bad # arguments");
639 greg 2.6 /* record surface struck */
640 greg 2.9 hitfront = (r->rod > 0);
641 greg 2.1 /* load cal file */
642 greg 2.50 mf = hasthick ? getfunc(m, 5, 0x1d, 1)
643     : getfunc(m, 4, 0xe, 1) ;
644 greg 2.25 setfunc(m, r);
645 greg 2.50 nd.thick = 0; /* set thickness */
646     if (hasthick) {
647     nd.thick = evalue(mf->ep[0]);
648     if ((-FTINY <= nd.thick) & (nd.thick <= FTINY))
649     nd.thick = 0;
650     }
651 greg 2.26 /* check backface visibility */
652     if (!hitfront & !backvis) {
653     raytrans(r);
654     return(1);
655     }
656 greg 2.5 /* check other rays to pass */
657 greg 2.34 if (nd.thick != 0 && (r->crtype & SHADOW ||
658     !(r->crtype & (SPECULAR|AMBIENT)) ||
659 greg 2.29 (nd.thick > 0) ^ hitfront)) {
660 greg 2.5 raytrans(r); /* hide our proxy */
661 greg 2.1 return(1);
662     }
663 greg 2.51 if (hasthick && r->crtype & SHADOW) /* early shadow check #1 */
664     return(1);
665 greg 2.31 nd.mp = m;
666     nd.pr = r;
667 greg 2.5 /* get BSDF data */
668 greg 2.50 nd.sd = loadBSDF(m->oargs.sarg[hasthick]);
669 greg 2.51 /* early shadow check #2 */
670 greg 2.55 if (r->crtype & SHADOW && (nd.sd->tf == NULL) & (nd.sd->tb == NULL)) {
671     SDfreeCache(nd.sd);
672 greg 2.34 return(1);
673 greg 2.55 }
674 greg 2.1 /* diffuse reflectance */
675 greg 2.6 if (hitfront) {
676 greg 2.31 cvt_sdcolor(nd.rdiff, &nd.sd->rLambFront);
677     if (m->oargs.nfargs >= 3) {
678     setcolor(ctmp, m->oargs.farg[0],
679 greg 2.1 m->oargs.farg[1],
680     m->oargs.farg[2]);
681 greg 2.31 addcolor(nd.rdiff, ctmp);
682     }
683 greg 2.1 } else {
684 greg 2.31 cvt_sdcolor(nd.rdiff, &nd.sd->rLambBack);
685     if (m->oargs.nfargs >= 6) {
686     setcolor(ctmp, m->oargs.farg[3],
687 greg 2.1 m->oargs.farg[4],
688     m->oargs.farg[5]);
689 greg 2.31 addcolor(nd.rdiff, ctmp);
690     }
691 greg 2.1 }
692     /* diffuse transmittance */
693 greg 2.31 cvt_sdcolor(nd.tdiff, &nd.sd->tLamb);
694     if (m->oargs.nfargs >= 9) {
695     setcolor(ctmp, m->oargs.farg[6],
696 greg 2.1 m->oargs.farg[7],
697     m->oargs.farg[8]);
698 greg 2.31 addcolor(nd.tdiff, ctmp);
699     }
700 greg 2.1 /* get modifiers */
701     raytexture(r, m->omod);
702     /* modify diffuse values */
703     multcolor(nd.rdiff, r->pcol);
704     multcolor(nd.tdiff, r->pcol);
705     /* get up vector */
706 greg 2.50 upvec[0] = evalue(mf->ep[hasthick+0]);
707     upvec[1] = evalue(mf->ep[hasthick+1]);
708     upvec[2] = evalue(mf->ep[hasthick+2]);
709 greg 2.1 /* return to world coords */
710 greg 2.21 if (mf->fxp != &unitxf) {
711     multv3(upvec, upvec, mf->fxp->xfm);
712     nd.thick *= mf->fxp->sca;
713 greg 2.1 }
714 greg 2.23 if (r->rox != NULL) {
715     multv3(upvec, upvec, r->rox->f.xfm);
716     nd.thick *= r->rox->f.sca;
717     }
718 greg 2.1 raynormal(nd.pnorm, r);
719     /* compute local BSDF xform */
720     ec = SDcompXform(nd.toloc, nd.pnorm, upvec);
721     if (!ec) {
722 greg 2.4 nd.vray[0] = -r->rdir[0];
723     nd.vray[1] = -r->rdir[1];
724     nd.vray[2] = -r->rdir[2];
725     ec = SDmapDir(nd.vray, nd.toloc, nd.vray);
726 greg 2.20 }
727 greg 2.19 if (ec) {
728     objerror(m, WARNING, "Illegal orientation vector");
729 greg 2.55 SDfreeCache(nd.sd);
730 greg 2.19 return(1);
731 greg 2.1 }
732 greg 2.50 setcolor(nd.cthru, 0, 0, 0); /* consider through component */
733 greg 2.52 if (m->otype == MAT_ABSDF) {
734 greg 2.50 compute_through(&nd);
735     if (r->crtype & SHADOW) {
736     RAY tr; /* attempt to pass shadow ray */
737 greg 2.55 SDfreeCache(nd.sd);
738 greg 2.50 if (rayorigin(&tr, TRANS, r, nd.cthru) < 0)
739     return(1); /* no through component */
740     VCOPY(tr.rdir, r->rdir);
741     rayvalue(&tr); /* transmit with scaling */
742     multcolor(tr.rcol, tr.rcoef);
743     copycolor(r->rcol, tr.rcol);
744     return(1); /* we're done */
745     }
746 greg 2.34 }
747     ec = SDinvXform(nd.fromloc, nd.toloc);
748     if (!ec) /* determine BSDF resolution */
749     ec = SDsizeBSDF(nd.sr_vpsa, nd.vray, NULL,
750     SDqueryMin+SDqueryMax, nd.sd);
751 greg 2.20 if (ec)
752     objerror(m, USER, transSDError(ec));
753    
754 greg 2.9 nd.sr_vpsa[0] = sqrt(nd.sr_vpsa[0]);
755     nd.sr_vpsa[1] = sqrt(nd.sr_vpsa[1]);
756 greg 2.6 if (!hitfront) { /* perturb normal towards hit */
757 greg 2.1 nd.pnorm[0] = -nd.pnorm[0];
758     nd.pnorm[1] = -nd.pnorm[1];
759     nd.pnorm[2] = -nd.pnorm[2];
760     }
761     /* sample reflection */
762     sample_sdf(&nd, SDsampSpR);
763     /* sample transmission */
764     sample_sdf(&nd, SDsampSpT);
765     /* compute indirect diffuse */
766 greg 2.39 copycolor(ctmp, nd.rdiff);
767     addcolor(ctmp, nd.runsamp);
768     if (bright(ctmp) > FTINY) { /* ambient from reflection */
769 greg 2.6 if (!hitfront)
770 greg 2.1 flipsurface(r);
771     multambient(ctmp, r, nd.pnorm);
772     addcolor(r->rcol, ctmp);
773 greg 2.6 if (!hitfront)
774 greg 2.1 flipsurface(r);
775     }
776 greg 2.39 copycolor(ctmp, nd.tdiff);
777     addcolor(ctmp, nd.tunsamp);
778     if (bright(ctmp) > FTINY) { /* ambient from other side */
779 greg 2.1 FVECT bnorm;
780 greg 2.6 if (hitfront)
781 greg 2.1 flipsurface(r);
782     bnorm[0] = -nd.pnorm[0];
783     bnorm[1] = -nd.pnorm[1];
784     bnorm[2] = -nd.pnorm[2];
785 greg 2.9 if (nd.thick != 0) { /* proxy with offset? */
786 greg 2.5 VCOPY(vtmp, r->rop);
787 greg 2.18 VSUM(r->rop, vtmp, r->ron, nd.thick);
788 greg 2.5 multambient(ctmp, r, bnorm);
789     VCOPY(r->rop, vtmp);
790     } else
791     multambient(ctmp, r, bnorm);
792 greg 2.1 addcolor(r->rcol, ctmp);
793 greg 2.6 if (hitfront)
794 greg 2.1 flipsurface(r);
795     }
796     /* add direct component */
797 greg 2.22 if ((bright(nd.tdiff) <= FTINY) & (nd.sd->tf == NULL) &
798     (nd.sd->tb == NULL)) {
799 greg 2.5 direct(r, dir_brdf, &nd); /* reflection only */
800 greg 2.9 } else if (nd.thick == 0) {
801 greg 2.5 direct(r, dir_bsdf, &nd); /* thin surface scattering */
802     } else {
803     direct(r, dir_brdf, &nd); /* reflection first */
804     VCOPY(vtmp, r->rop); /* offset for transmitted */
805     VSUM(r->rop, vtmp, r->ron, -nd.thick);
806 greg 2.6 direct(r, dir_btdf, &nd); /* separate transmission */
807 greg 2.5 VCOPY(r->rop, vtmp);
808     }
809 greg 2.1 /* clean up */
810     SDfreeCache(nd.sd);
811     return(1);
812     }