--- ray/src/rt/m_bsdf.c 2011/02/18 00:40:25 2.1 +++ ray/src/rt/m_bsdf.c 2012/04/24 15:36:40 2.19 @@ -1,5 +1,5 @@ #ifndef lint -static const char RCSid[] = "$Id: m_bsdf.c,v 2.1 2011/02/18 00:40:25 greg Exp $"; +static const char RCSid[] = "$Id: m_bsdf.c,v 2.19 2012/04/24 15:36:40 greg Exp $"; #endif /* * Shading for materials with BSDFs taken from XML data files @@ -8,7 +8,6 @@ static const char RCSid[] = "$Id: m_bsdf.c,v 2.1 2011/ #include "copyright.h" #include "ray.h" -#include "paths.h" #include "ambient.h" #include "source.h" #include "func.h" @@ -18,41 +17,57 @@ static const char RCSid[] = "$Id: m_bsdf.c,v 2.1 2011/ /* * Arguments to this material include optional diffuse colors. * String arguments include the BSDF and function files. - * A thickness variable causes the strange but useful behavior - * of translating transmitted rays this distance past the surface - * intersection in the normal direction to bypass intervening geometry. - * This only affects scattered, non-source directed samples. Thus, - * thickness is relevant only if there is a transmitted component. - * A positive thickness has the further side-effect that an unscattered + * A non-zero thickness causes the strange but useful behavior + * of translating transmitted rays this distance beneath the surface + * (opposite the surface normal) to bypass any intervening geometry. + * Translation only affects scattered, non-source-directed samples. + * A non-zero thickness has the further side-effect that an unscattered * (view) ray will pass right through our material if it has any - * non-diffuse transmission, making our BSDF invisible. This allows the - * underlying geometry to become visible. A matching surface should be - * placed on the other side, less than the thickness away, if the backside - * reflectance is non-zero. + * non-diffuse transmission, making the BSDF surface invisible. This + * shows the proxied geometry instead. Thickness has the further + * effect of turning off reflection on the hidden side so that rays + * heading in the opposite direction pass unimpeded through the BSDF + * surface. A paired surface may be placed on the opposide side of + * the detail geometry, less than this thickness away, if a two-way + * proxy is desired. Note that the sign of the thickness is important. + * A positive thickness hides geometry behind the BSDF surface and uses + * front reflectance and transmission properties. A negative thickness + * hides geometry in front of the surface when rays hit from behind, + * and applies only the transmission and backside reflectance properties. + * Reflection is ignored on the hidden side, as those rays pass through. * The "up" vector for the BSDF is given by three variables, defined * (along with the thickness) by the named function file, or '.' if none. * Together with the surface normal, this defines the local coordinate * system for the BSDF. * We do not reorient the surface, so if the BSDF has no back-side - * reflectance and none is given in the real arguments, the surface will - * appear as black when viewed from behind (unless backvis is false). - * The diffuse compnent arguments are added to components in the BSDF file, + * reflectance and none is given in the real arguments, a BSDF surface + * with zero thickness will appear black when viewed from behind + * unless backface visibility is off. + * The diffuse arguments are added to components in the BSDF file, * not multiplied. However, patterns affect this material as a multiplier * on everything except non-diffuse reflection. * * Arguments for MAT_BSDF are: * 6+ thick BSDFfile ux uy uz funcfile transform * 0 - * 0|3|9 rdf gdf bdf + * 0|3|6|9 rdf gdf bdf * rdb gdb bdb * rdt gdt bdt */ +/* + * Note that our reverse ray-tracing process means that the positions + * of incoming and outgoing vectors may be reversed in our calls + * to the BSDF library. This is fine, since the bidirectional nature + * of the BSDF (that's what the 'B' stands for) means it all works out. + */ + typedef struct { OBJREC *mp; /* material pointer */ RAY *pr; /* intersected ray */ FVECT pnorm; /* perturbed surface normal */ - FVECT vinc; /* local incident vector */ + FVECT vray; /* local outgoing (return) vector */ + double sr_vpsa[2]; /* sqrt of BSDF projected solid angle extrema */ RREAL toloc[3][3]; /* world to local BSDF coords */ RREAL fromloc[3][3]; /* local BSDF coords to world */ double thick; /* surface thickness */ @@ -65,29 +80,96 @@ typedef struct { #define cvt_sdcolor(cv, svp) ccy2rgb(&(svp)->spec, (svp)->cieY, cv) -/* Convert error from BSDF library */ -static char * -cvt_sderr(SDError ec) +/* Jitter ray sample according to projected solid angle and specjitter */ +static void +bsdf_jitter(FVECT vres, BSDFDAT *ndp, double sr_psa) { - if (!SDerrorDetail[0]) - return(strcpy(errmsg, SDerrorEnglish[ec])); - sprintf(errmsg, "%s: %s", SDerrorEnglish[ec], SDerrorDetail); - return(errmsg); + VCOPY(vres, ndp->vray); + if (specjitter < 1.) + sr_psa *= specjitter; + if (sr_psa <= FTINY) + return; + vres[0] += sr_psa*(.5 - frandom()); + vres[1] += sr_psa*(.5 - frandom()); + normalize(vres); } -/* Compute source contribution for BSDF */ +/* Evaluate BSDF for direct component, returning true if OK to proceed */ +static int +direct_bsdf_OK(COLOR cval, FVECT ldir, double omega, BSDFDAT *ndp) +{ + int nsamp, ok = 0; + FVECT vsrc, vsmp, vjit; + double tomega; + double sf, tsr, sd[2]; + COLOR csmp; + SDValue sv; + SDError ec; + int i; + /* transform source direction */ + if (SDmapDir(vsrc, ndp->toloc, ldir) != SDEnone) + return(0); + /* assign number of samples */ + ec = SDsizeBSDF(&tomega, ndp->vray, vsrc, SDqueryMin, ndp->sd); + if (ec) + goto baderror; + /* check indirect over-counting */ + if (ndp->thick != 0 && ndp->pr->crtype & (SPECULAR|AMBIENT) + && vsrc[2] > 0 ^ ndp->vray[2] > 0) { + double dx = vsrc[0] + ndp->vray[0]; + double dy = vsrc[1] + ndp->vray[1]; + if (dx*dx + dy*dy <= omega+tomega) + return(0); + } + sf = specjitter * ndp->pr->rweight; + if (25.*tomega <= omega) + nsamp = 100.*sf + .5; + else + nsamp = 4.*sf*omega/tomega + .5; + nsamp += !nsamp; + setcolor(cval, .0, .0, .0); /* sample our source area */ + sf = sqrt(omega); + tsr = sqrt(tomega); + for (i = nsamp; i--; ) { + VCOPY(vsmp, vsrc); /* jitter query directions */ + if (nsamp > 1) { + multisamp(sd, 2, (i + frandom())/(double)nsamp); + vsmp[0] += (sd[0] - .5)*sf; + vsmp[1] += (sd[1] - .5)*sf; + if (normalize(vsmp) == 0) { + --nsamp; + continue; + } + } + bsdf_jitter(vjit, ndp, tsr); + /* compute BSDF */ + ec = SDevalBSDF(&sv, vjit, vsmp, ndp->sd); + if (ec) + goto baderror; + if (sv.cieY <= FTINY) /* worth using? */ + continue; + cvt_sdcolor(csmp, &sv); + addcolor(cval, csmp); /* average it in */ + ++ok; + } + sf = 1./(double)nsamp; + scalecolor(cval, sf); + return(ok); +baderror: + objerror(ndp->mp, USER, transSDError(ec)); + return(0); /* gratis return */ +} + +/* Compute source contribution for BSDF (reflected & transmitted) */ static void -dirbsdf( +dir_bsdf( COLOR cval, /* returned coefficient */ void *nnp, /* material data */ FVECT ldir, /* light source direction */ double omega /* light source size */ ) { - BSDFDAT *np = nnp; - SDError ec; - SDValue sv; - FVECT vout; + BSDFDAT *np = (BSDFDAT *)nnp; double ldot; double dtmp; COLOR ctmp; @@ -98,7 +180,7 @@ dirbsdf( if ((-FTINY <= ldot) & (ldot <= FTINY)) return; - if (ldot > .0 && bright(np->rdiff) > FTINY) { + if (ldot > 0 && bright(np->rdiff) > FTINY) { /* * Compute added diffuse reflected component. */ @@ -107,7 +189,7 @@ dirbsdf( scalecolor(ctmp, dtmp); addcolor(cval, ctmp); } - if (ldot < .0 && bright(np->tdiff) > FTINY) { + if (ldot < 0 && bright(np->tdiff) > FTINY) { /* * Compute added diffuse transmission. */ @@ -119,25 +201,19 @@ dirbsdf( /* * Compute scattering coefficient using BSDF. */ - if (SDmapDir(vout, np->toloc, ldir) != SDEnone) + if (!direct_bsdf_OK(ctmp, ldir, omega, np)) return; - ec = SDevalBSDF(&sv, vout, np->vinc, np->sd); - if (ec) - objerror(np->mp, USER, cvt_sderr(ec)); - - if (sv.cieY <= FTINY) /* not worth using? */ - return; - cvt_sdcolor(ctmp, &sv); - if (ldot > .0) { /* pattern only diffuse reflection */ + if (ldot > 0) { /* pattern only diffuse reflection */ COLOR ctmp1, ctmp2; - dtmp = (np->pr->rod > .0) ? np->sd->rLambFront.cieY + dtmp = (np->pr->rod > 0) ? np->sd->rLambFront.cieY : np->sd->rLambBack.cieY; - dtmp /= PI * sv.cieY; /* diffuse fraction */ + /* diffuse fraction */ + dtmp /= PI * bright(ctmp); copycolor(ctmp2, np->pr->pcol); scalecolor(ctmp2, dtmp); setcolor(ctmp1, 1.-dtmp, 1.-dtmp, 1.-dtmp); addcolor(ctmp1, ctmp2); - multcolor(ctmp, ctmp1); /* apply desaturated pattern */ + multcolor(ctmp, ctmp1); /* apply derated pattern */ dtmp = ldot * omega; } else { /* full pattern on transmission */ multcolor(ctmp, np->pr->pcol); @@ -147,64 +223,149 @@ dirbsdf( addcolor(cval, ctmp); } +/* Compute source contribution for BSDF (reflected only) */ +static void +dir_brdf( + COLOR cval, /* returned coefficient */ + void *nnp, /* material data */ + FVECT ldir, /* light source direction */ + double omega /* light source size */ +) +{ + BSDFDAT *np = (BSDFDAT *)nnp; + double ldot; + double dtmp; + COLOR ctmp, ctmp1, ctmp2; + + setcolor(cval, .0, .0, .0); + + ldot = DOT(np->pnorm, ldir); + + if (ldot <= FTINY) + return; + + if (bright(np->rdiff) > FTINY) { + /* + * Compute added diffuse reflected component. + */ + copycolor(ctmp, np->rdiff); + dtmp = ldot * omega * (1./PI); + scalecolor(ctmp, dtmp); + addcolor(cval, ctmp); + } + /* + * Compute reflection coefficient using BSDF. + */ + if (!direct_bsdf_OK(ctmp, ldir, omega, np)) + return; + /* pattern only diffuse reflection */ + dtmp = (np->pr->rod > 0) ? np->sd->rLambFront.cieY + : np->sd->rLambBack.cieY; + dtmp /= PI * bright(ctmp); /* diffuse fraction */ + copycolor(ctmp2, np->pr->pcol); + scalecolor(ctmp2, dtmp); + setcolor(ctmp1, 1.-dtmp, 1.-dtmp, 1.-dtmp); + addcolor(ctmp1, ctmp2); + multcolor(ctmp, ctmp1); /* apply derated pattern */ + dtmp = ldot * omega; + scalecolor(ctmp, dtmp); + addcolor(cval, ctmp); +} + +/* Compute source contribution for BSDF (transmitted only) */ +static void +dir_btdf( + COLOR cval, /* returned coefficient */ + void *nnp, /* material data */ + FVECT ldir, /* light source direction */ + double omega /* light source size */ +) +{ + BSDFDAT *np = (BSDFDAT *)nnp; + double ldot; + double dtmp; + COLOR ctmp; + + setcolor(cval, .0, .0, .0); + + ldot = DOT(np->pnorm, ldir); + + if (ldot >= -FTINY) + return; + + if (bright(np->tdiff) > FTINY) { + /* + * Compute added diffuse transmission. + */ + copycolor(ctmp, np->tdiff); + dtmp = -ldot * omega * (1.0/PI); + scalecolor(ctmp, dtmp); + addcolor(cval, ctmp); + } + /* + * Compute scattering coefficient using BSDF. + */ + if (!direct_bsdf_OK(ctmp, ldir, omega, np)) + return; + /* full pattern on transmission */ + multcolor(ctmp, np->pr->pcol); + dtmp = -ldot * omega; + scalecolor(ctmp, dtmp); + addcolor(cval, ctmp); +} + /* Sample separate BSDF component */ static int sample_sdcomp(BSDFDAT *ndp, SDComponent *dcp, int usepat) { int nstarget = 1; - int nsent = 0; + int nsent; SDError ec; SDValue bsv; - double sthick; - FVECT vout; + double xrand; + FVECT vsmp; RAY sr; - int ntrials; /* multiple samples? */ if (specjitter > 1.5) { nstarget = specjitter*ndp->pr->rweight + .5; - if (nstarget < 1) - nstarget = 1; + nstarget += !nstarget; } - /* run through our trials */ - for (ntrials = 0; nsent < nstarget && ntrials < 9*nstarget; ntrials++) { - SDerrorDetail[0] = '\0'; - /* sample direction & coef. */ - ec = SDsampComponent(&bsv, vout, ndp->vinc, - ntrials ? frandom() - : urand(ilhash(dimlist,ndims)+samplendx), - dcp); + /* run through our samples */ + for (nsent = 0; nsent < nstarget; nsent++) { + if (nstarget == 1) { /* stratify random variable */ + xrand = urand(ilhash(dimlist,ndims)+samplendx); + if (specjitter < 1.) + xrand = .5 + specjitter*(xrand-.5); + } else { + xrand = (nsent + frandom())/(double)nstarget; + } + SDerrorDetail[0] = '\0'; /* sample direction & coef. */ + bsdf_jitter(vsmp, ndp, ndp->sr_vpsa[0]); + ec = SDsampComponent(&bsv, vsmp, xrand, dcp); if (ec) - objerror(ndp->mp, USER, cvt_sderr(ec)); - /* zero component? */ - if (bsv.cieY <= FTINY) + objerror(ndp->mp, USER, transSDError(ec)); + if (bsv.cieY <= FTINY) /* zero component? */ break; /* map vector to world */ - if (SDmapDir(sr.rdir, ndp->fromloc, vout) != SDEnone) + if (SDmapDir(sr.rdir, ndp->fromloc, vsmp) != SDEnone) break; - /* unintentional penetration? */ - if (DOT(sr.rdir, ndp->pr->ron) > .0 ^ vout[2] > .0) - continue; /* spawn a specular ray */ if (nstarget > 1) bsv.cieY /= (double)nstarget; - cvt_sdcolor(sr.rcoef, &bsv); /* use color */ - if (usepat) /* pattern on transmission */ + cvt_sdcolor(sr.rcoef, &bsv); /* use sample color */ + if (usepat) /* apply pattern? */ multcolor(sr.rcoef, ndp->pr->pcol); if (rayorigin(&sr, SPECULAR, ndp->pr, sr.rcoef) < 0) { - if (maxdepth > 0) + if (maxdepth > 0) break; - ++nsent; /* Russian roulette victim */ - continue; + continue; /* Russian roulette victim */ } - /* need to move origin? */ - sthick = (ndp->pr->rod > .0) ? -ndp->thick : ndp->thick; - if (sthick < .0 ^ vout[2] > .0) - VSUM(sr.rorg, sr.rorg, ndp->pr->ron, sthick); - + /* need to offset origin? */ + if (ndp->thick != 0 && ndp->pr->rod > 0 ^ vsmp[2] > 0) + VSUM(sr.rorg, sr.rorg, ndp->pr->ron, -ndp->thick); rayvalue(&sr); /* send & evaluate sample */ multcolor(sr.rcol, sr.rcoef); addcolor(ndp->pr->rcol, sr.rcol); - ++nsent; } return(nsent); } @@ -223,7 +384,7 @@ sample_sdf(BSDFDAT *ndp, int sflags) cvt_sdcolor(unsc, &ndp->sd->tLamb); } else /* sflags == SDsampSpR */ { unsc = ndp->runsamp; - if (ndp->pr->rod > .0) { + if (ndp->pr->rod > 0) { dfp = ndp->sd->rf; cvt_sdcolor(unsc, &ndp->sd->rLambFront); } else { @@ -236,9 +397,12 @@ sample_sdf(BSDFDAT *ndp, int sflags) return(0); /* below sampling threshold? */ if (dfp->maxHemi <= specthresh+FTINY) { - if (dfp->maxHemi > FTINY) { /* XXX no color from BSDF! */ - double d = SDdirectHemi(ndp->vinc, sflags, ndp->sd); + if (dfp->maxHemi > FTINY) { /* XXX no color from BSDF */ + FVECT vjit; + double d; COLOR ctmp; + bsdf_jitter(vjit, ndp, ndp->sr_vpsa[1]); + d = SDdirectHemi(vjit, sflags, ndp->sd); if (sflags == SDsampSpT) { copycolor(ctmp, ndp->pr->pcol); scalecolor(ctmp, d); @@ -263,55 +427,40 @@ sample_sdf(BSDFDAT *ndp, int sflags) int m_bsdf(OBJREC *m, RAY *r) { + int hitfront; COLOR ctmp; SDError ec; - FVECT upvec, outVec; + FVECT upvec, vtmp; MFUNC *mf; BSDFDAT nd; /* check arguments */ if ((m->oargs.nsargs < 6) | (m->oargs.nfargs > 9) | (m->oargs.nfargs % 3)) objerror(m, USER, "bad # arguments"); - - SDerrorDetail[0] = '\0'; /* get BSDF data */ - nd.sd = SDgetCache(m->oargs.sarg[1]); - if (nd.sd == NULL) - error(SYSTEM, "out of memory in m_bsdf"); - if (!SDisLoaded(nd.sd)) { - char *pname = getpath(m->oargs.sarg[1], getrlibpath(), R_OK); - if (pname == NULL) { - sprintf(errmsg, "cannot find BSDF file \"%s\"", - m->oargs.sarg[1]); - objerror(m, USER, errmsg); - } - ec = SDloadFile(nd.sd, pname); - if (ec) - objerror(m, USER, cvt_sderr(ec)); - SDretainSet = SDretainAll; - } + /* record surface struck */ + hitfront = (r->rod > 0); /* load cal file */ mf = getfunc(m, 5, 0x1d, 1); /* get thickness */ nd.thick = evalue(mf->ep[0]); - if (nd.thick < .0) + if ((-FTINY <= nd.thick) & (nd.thick <= FTINY)) nd.thick = .0; /* check shadow */ if (r->crtype & SHADOW) { - SDfreeCache(nd.sd); - if (nd.thick > FTINY && nd.sd->tf != NULL && - nd.sd->tf->maxHemi > FTINY) + if (nd.thick != 0) raytrans(r); /* pass-through */ - return(1); /* else shadow */ + return(1); /* or shadow */ } - /* check unscattered ray */ - if (!(r->crtype & (SPECULAR|AMBIENT)) && nd.thick > FTINY && - nd.sd->tf != NULL && nd.sd->tf->maxHemi > FTINY) { - SDfreeCache(nd.sd); - raytrans(r); /* pass-through */ + /* check other rays to pass */ + if (nd.thick != 0 && (!(r->crtype & (SPECULAR|AMBIENT)) || + nd.thick > 0 ^ hitfront)) { + raytrans(r); /* hide our proxy */ return(1); } + /* get BSDF data */ + nd.sd = loadBSDF(m->oargs.sarg[1]); /* diffuse reflectance */ - if (r->rod > .0) { + if (hitfront) { if (m->oargs.nfargs < 3) setcolor(nd.rdiff, .0, .0, .0); else @@ -320,10 +469,8 @@ m_bsdf(OBJREC *m, RAY *r) m->oargs.farg[2]); } else { if (m->oargs.nfargs < 6) { /* check invisible backside */ - if (!backvis && (nd.sd->rb == NULL || - nd.sd->rb->maxHemi <= FTINY) && - (nd.sd->tf == NULL || - nd.sd->tf->maxHemi <= FTINY)) { + if (!backvis && (nd.sd->rb == NULL) & + (nd.sd->tf == NULL)) { SDfreeCache(nd.sd); raytrans(r); return(1); @@ -345,10 +492,6 @@ m_bsdf(OBJREC *m, RAY *r) nd.pr = r; /* get modifiers */ raytexture(r, m->omod); - if (bright(r->pcol) <= FTINY) { /* black pattern?! */ - SDfreeCache(nd.sd); - return(1); - } /* modify diffuse values */ multcolor(nd.rdiff, r->pcol); multcolor(nd.tdiff, r->pcol); @@ -365,19 +508,28 @@ m_bsdf(OBJREC *m, RAY *r) /* compute local BSDF xform */ ec = SDcompXform(nd.toloc, nd.pnorm, upvec); if (!ec) { - nd.vinc[0] = -r->rdir[0]; - nd.vinc[1] = -r->rdir[1]; - nd.vinc[2] = -r->rdir[2]; - ec = SDmapDir(nd.vinc, nd.toloc, nd.vinc); + nd.vray[0] = -r->rdir[0]; + nd.vray[1] = -r->rdir[1]; + nd.vray[2] = -r->rdir[2]; + ec = SDmapDir(nd.vray, nd.toloc, nd.vray); + } + if (ec) { + objerror(m, WARNING, "Illegal orientation vector"); + return(1); } + ec = SDinvXform(nd.fromloc, nd.toloc); + /* determine BSDF resolution */ if (!ec) - ec = SDinvXform(nd.fromloc, nd.toloc); + ec = SDsizeBSDF(nd.sr_vpsa, nd.vray, NULL, + SDqueryMin+SDqueryMax, nd.sd); if (ec) { - objerror(m, WARNING, cvt_sderr(ec)); + objerror(m, WARNING, transSDError(ec)); SDfreeCache(nd.sd); return(1); } - if (r->rod < .0) { /* perturb normal towards hit */ + nd.sr_vpsa[0] = sqrt(nd.sr_vpsa[0]); + nd.sr_vpsa[1] = sqrt(nd.sr_vpsa[1]); + if (!hitfront) { /* perturb normal towards hit */ nd.pnorm[0] = -nd.pnorm[0]; nd.pnorm[1] = -nd.pnorm[1]; nd.pnorm[2] = -nd.pnorm[2]; @@ -389,30 +541,46 @@ m_bsdf(OBJREC *m, RAY *r) /* compute indirect diffuse */ copycolor(ctmp, nd.rdiff); addcolor(ctmp, nd.runsamp); - if (bright(ctmp) > FTINY) { /* ambient from this side */ - if (r->rod < .0) + if (bright(ctmp) > FTINY) { /* ambient from reflection */ + if (!hitfront) flipsurface(r); multambient(ctmp, r, nd.pnorm); addcolor(r->rcol, ctmp); - if (r->rod < .0) + if (!hitfront) flipsurface(r); } copycolor(ctmp, nd.tdiff); addcolor(ctmp, nd.tunsamp); if (bright(ctmp) > FTINY) { /* ambient from other side */ FVECT bnorm; - if (r->rod > .0) + if (hitfront) flipsurface(r); bnorm[0] = -nd.pnorm[0]; bnorm[1] = -nd.pnorm[1]; bnorm[2] = -nd.pnorm[2]; - multambient(ctmp, r, bnorm); + if (nd.thick != 0) { /* proxy with offset? */ + VCOPY(vtmp, r->rop); + VSUM(r->rop, vtmp, r->ron, nd.thick); + multambient(ctmp, r, bnorm); + VCOPY(r->rop, vtmp); + } else + multambient(ctmp, r, bnorm); addcolor(r->rcol, ctmp); - if (r->rod > .0) + if (hitfront) flipsurface(r); } /* add direct component */ - direct(r, dirbsdf, &nd); + if ((bright(nd.tdiff) <= FTINY) & (nd.sd->tf == NULL)) { + direct(r, dir_brdf, &nd); /* reflection only */ + } else if (nd.thick == 0) { + direct(r, dir_bsdf, &nd); /* thin surface scattering */ + } else { + direct(r, dir_brdf, &nd); /* reflection first */ + VCOPY(vtmp, r->rop); /* offset for transmitted */ + VSUM(r->rop, vtmp, r->ron, -nd.thick); + direct(r, dir_btdf, &nd); /* separate transmission */ + VCOPY(r->rop, vtmp); + } /* clean up */ SDfreeCache(nd.sd); return(1);