--- ray/src/rt/dielectric.c 1989/06/07 08:38:31 1.2 +++ ray/src/rt/dielectric.c 1996/04/17 14:01:52 2.11 @@ -16,6 +16,8 @@ static char SCCSid[] = "$SunId$ LBL"; #ifdef DISPERSE #include "source.h" +static disperse(); +static int lambda(); #endif /* @@ -48,14 +50,29 @@ static char SCCSid[] = "$SunId$ LBL"; #define MINCOS 0.997 /* minimum dot product for dispersion */ +extern COLOR cextinction; /* global coefficient of extinction */ +extern COLOR salbedo; /* global scattering albedo */ -m_dielectric(m, r) /* color a ray which hit something transparent */ + +static double +mylog(x) /* special log for extinction coefficients */ +double x; +{ + if (x < 1e-40) + return(-100.); + if (x >= 1.) + return(0.); + return(log(x)); +} + + +m_dielectric(m, r) /* color a ray which hit a dielectric interface */ OBJREC *m; register RAY *r; { - double sqrt(), pow(); double cos1, cos2, nratio; - COLOR mcolor; + COLOR ctrans; + COLOR talb; double mabsorp; double refl, trans; FVECT dnorm; @@ -80,19 +97,38 @@ register RAY *r; dnorm[0] = -dnorm[0]; dnorm[1] = -dnorm[1]; dnorm[2] = -dnorm[2]; - setcolor(mcolor, pow(m->oargs.farg[0], r->rot), - pow(m->oargs.farg[1], r->rot), - pow(m->oargs.farg[2], r->rot)); + setcolor(r->cext, -mylog(m->oargs.farg[0]*colval(r->pcol,RED)), + -mylog(m->oargs.farg[1]*colval(r->pcol,GRN)), + -mylog(m->oargs.farg[2]*colval(r->pcol,BLU))); + setcolor(r->albedo, 0., 0., 0.); + r->gecc = 0.; + if (m->otype == MAT_INTERFACE) { + setcolor(ctrans, + -mylog(m->oargs.farg[4]*colval(r->pcol,RED)), + -mylog(m->oargs.farg[5]*colval(r->pcol,GRN)), + -mylog(m->oargs.farg[6]*colval(r->pcol,BLU))); + setcolor(talb, 0., 0., 0.); + } else { + copycolor(ctrans, cextinction); + copycolor(talb, salbedo); + } } else { /* outside */ nratio = 1.0 / nratio; - if (m->otype == MAT_INTERFACE) - setcolor(mcolor, pow(m->oargs.farg[4], r->rot), - pow(m->oargs.farg[5], r->rot), - pow(m->oargs.farg[6], r->rot)); - else - setcolor(mcolor, 1.0, 1.0, 1.0); + + setcolor(ctrans, -mylog(m->oargs.farg[0]*colval(r->pcol,RED)), + -mylog(m->oargs.farg[1]*colval(r->pcol,GRN)), + -mylog(m->oargs.farg[2]*colval(r->pcol,BLU))); + setcolor(talb, 0., 0., 0.); + if (m->otype == MAT_INTERFACE) { + setcolor(r->cext, + -mylog(m->oargs.farg[4]*colval(r->pcol,RED)), + -mylog(m->oargs.farg[5]*colval(r->pcol,GRN)), + -mylog(m->oargs.farg[6]*colval(r->pcol,BLU))); + setcolor(r->albedo, 0., 0., 0.); + r->gecc = 0.; + } } - mabsorp = bright(mcolor); + mabsorp = exp(-bright(r->cext)*r->rot); /* approximate */ d2 = 1.0 - nratio*nratio*(1.0 - cos1*cos1); /* compute cos theta2 */ @@ -113,7 +149,7 @@ register RAY *r; d1 = (d1 - d2) / (d1 + d2); refl += d1 * d1; - refl /= 2.0; + refl *= 0.5; trans = 1.0 - refl; if (rayorigin(&p, r, REFRACTED, mabsorp*trans) == 0) { @@ -127,14 +163,18 @@ register RAY *r; if (m->otype != MAT_DIELECTRIC || r->rod > 0.0 || r->crtype & SHADOW + || !directvis || m->oargs.farg[4] == 0.0 || !disperse(m, r, p.rdir, trans)) #endif { + copycolor(p.cext, ctrans); + copycolor(p.albedo, talb); rayvalue(&p); - multcolor(mcolor, r->pcol); /* modify */ scalecolor(p.rcol, trans); addcolor(r->rcol, p.rcol); + if (nratio >= 1.0-FTINY && nratio <= 1.0+FTINY) + r->rt = r->rot + p.rt; } } } @@ -151,8 +191,8 @@ register RAY *r; scalecolor(p.rcol, refl); /* color contribution */ addcolor(r->rcol, p.rcol); } - - multcolor(r->rcol, mcolor); /* multiply by transmittance */ + /* rayvalue() computes absorption */ + return(1); } @@ -165,13 +205,12 @@ RAY *r; FVECT vt; double tr; { - double sqrt(); RAY sray, *entray; FVECT v1, v2, n1, n2; FVECT dv, v2Xdv; double v2Xdvv2Xdv; - int sn, success = 0; - double omega; + int success = 0; + SRCINDEX si; FVECT vtmp1, vtmp2; double dtmp1, dtmp2; int l1, l2; @@ -233,12 +272,11 @@ double tr; v2Xdvv2Xdv = DOT(v2Xdv, v2Xdv); /* check sources */ - for (sn = 0; sn < nsources; sn++) { + initsrcindex(&si); + while (srcray(&sray, r, &si)) { - if ((omega = srcray(&sray, r, sn)) == 0.0 || - DOT(sray.rdir, v2) < MINCOS) + if (DOT(sray.rdir, v2) < MINCOS) continue; /* bad source */ - /* adjust source ray */ dtmp1 = DOT(v2Xdv, sray.rdir) / v2Xdvv2Xdv; @@ -263,7 +301,7 @@ double tr; */ fcross(vtmp1, v2Xdv, sray.rdir); - dtmp1 = sqrt(omega / v2Xdvv2Xdv / PI); + dtmp1 = sqrt(si.dom / v2Xdvv2Xdv / PI); /* compute first ray */ for (i = 0; i < 3; i++)