--- ray/src/common/bsdf.c	2011/04/17 17:45:13	2.21
+++ ray/src/common/bsdf.c	2011/04/24 20:16:52	2.25
@@ -1,5 +1,5 @@
 #ifndef lint
-static const char RCSid[] = "$Id: bsdf.c,v 2.21 2011/04/17 17:45:13 greg Exp $";
+static const char RCSid[] = "$Id: bsdf.c,v 2.25 2011/04/24 20:16:52 greg Exp $";
 #endif
 /*
  *  bsdf.c
@@ -99,7 +99,7 @@ SDloadGeometry(SDData *sd, ezxml_t wdb)
 	if ((geom = ezxml_child(wdb, "Thickness")) != NULL)
 		sd->dim[2] = atof(ezxml_txt(geom)) *
 				to_meters(ezxml_attr(geom, "unit"));
-	if ((sd->dim[0] < .0) | (sd->dim[1] < .0) | (sd->dim[2] < .0)) {
+	if ((sd->dim[0] < 0) | (sd->dim[1] < 0) | (sd->dim[2] < 0)) {
 		sprintf(SDerrorDetail, "Negative size in \"%s\"", sd->name);
 		return SDEdata;
 	}
@@ -405,34 +405,35 @@ SDfreeCache(const SDData *sd)
 
 /* Sample an individual BSDF component */
 SDError
-SDsampComponent(SDValue *sv, FVECT outVec, const FVECT inVec,
-			double randX, SDComponent *sdc)
+SDsampComponent(SDValue *sv, FVECT ioVec, double randX, SDComponent *sdc)
 {
 	float		coef[SDmaxCh];
 	SDError		ec;
+	FVECT		inVec;
 	const SDCDst	*cd;
 	double		d;
 	int		n;
 					/* check arguments */
-	if ((sv == NULL) | (outVec == NULL) | (inVec == NULL) | (sdc == NULL))
+	if ((sv == NULL) | (ioVec == NULL) | (sdc == NULL))
 		return SDEargument;
 					/* get cumulative distribution */
+	VCOPY(inVec, ioVec);
 	cd = (*sdc->func->getCDist)(inVec, sdc);
 	if (cd == NULL)
 		return SDEmemory;
 	if (cd->cTotal <= 1e-7) {	/* anything to sample? */
 		sv->spec = c_dfcolor;
 		sv->cieY = .0;
-		memset(outVec, 0, 3*sizeof(double));
+		memset(ioVec, 0, 3*sizeof(double));
 		return SDEnone;
 	}
 	sv->cieY = cd->cTotal;
 					/* compute sample direction */
-	ec = (*sdc->func->sampCDist)(outVec, randX, cd);
+	ec = (*sdc->func->sampCDist)(ioVec, randX, cd);
 	if (ec)
 		return ec;
 					/* get BSDF color */
-	n = (*sdc->func->getBSDFs)(coef, outVec, inVec, sdc->dist);
+	n = (*sdc->func->getBSDFs)(coef, ioVec, inVec, sdc);
 	if (n <= 0) {
 		strcpy(SDerrorDetail, "BSDF sample value error");
 		return SDEinternal;
@@ -480,7 +481,7 @@ SDdiffuseSamp(FVECT outVec, int outFront, double randX
 	SDmultiSamp(outVec, 2, randX);
 	SDsquare2disk(outVec, outVec[0], outVec[1]);
 	outVec[2] = 1. - outVec[0]*outVec[0] - outVec[1]*outVec[1];
-	if (outVec[2] > .0)		/* a bit of paranoia */
+	if (outVec[2] > 0)		/* a bit of paranoia */
 		outVec[2] = sqrt(outVec[2]);
 	if (!outFront)			/* going out back? */
 		outVec[2] = -outVec[2];
@@ -488,13 +489,14 @@ SDdiffuseSamp(FVECT outVec, int outFront, double randX
 
 /* Query projected solid angle coverage for non-diffuse BSDF direction */
 SDError
-SDsizeBSDF(double *projSA, const FVECT vec, int qflags, const SDData *sd)
+SDsizeBSDF(double *projSA, const FVECT v1, const RREAL *v2,
+				int qflags, const SDData *sd)
 {
-	SDSpectralDF	*rdf;
+	SDSpectralDF	*rdf, *tdf;
 	SDError		ec;
 	int		i;
 					/* check arguments */
-	if ((projSA == NULL) | (vec == NULL) | (sd == NULL))
+	if ((projSA == NULL) | (v1 == NULL))
 		return SDEargument;
 					/* initialize extrema */
 	switch (qflags) {
@@ -510,20 +512,25 @@ SDsizeBSDF(double *projSA, const FVECT vec, int qflags
 	case 0:
 		return SDEargument;
 	}
-	if (vec[2] > .0)		/* front surface query? */
+	if (v1[2] > 0)			/* front surface query? */
 		rdf = sd->rf;
 	else
 		rdf = sd->rb;
+	tdf = NULL;			/* transmitted component? */
+	if (v2 != NULL && v1[2] > 0 ^ v2[2] > 0) {
+		rdf = NULL;
+		tdf = sd->tf;
+	}
 	ec = SDEdata;			/* run through components */
 	for (i = (rdf==NULL) ? 0 : rdf->ncomp; i--; ) {
-		ec = (*rdf->comp[i].func->queryProjSA)(projSA, vec, qflags,
-							rdf->comp[i].dist);
+		ec = (*rdf->comp[i].func->queryProjSA)(projSA, v1, v2,
+						qflags, &rdf->comp[i]);
 		if (ec)
 			return ec;
 	}
-	for (i = (sd->tf==NULL) ? 0 : sd->tf->ncomp; i--; ) {
-		ec = (*sd->tf->comp[i].func->queryProjSA)(projSA, vec, qflags,
-							sd->tf->comp[i].dist);
+	for (i = (tdf==NULL) ? 0 : tdf->ncomp; i--; ) {
+		ec = (*tdf->comp[i].func->queryProjSA)(projSA, v1, v2,
+						qflags, &tdf->comp[i]);
 		if (ec)
 			return ec;
 	}
@@ -547,8 +554,8 @@ SDevalBSDF(SDValue *sv, const FVECT outVec, const FVEC
 	if ((sv == NULL) | (outVec == NULL) | (inVec == NULL) | (sd == NULL))
 		return SDEargument;
 					/* whose side are we on? */
-	inFront = (inVec[2] > .0);
-	outFront = (outVec[2] > .0);
+	inFront = (inVec[2] > 0);
+	outFront = (outVec[2] > 0);
 					/* start with diffuse portion */
 	if (inFront & outFront) {
 		*sv = sd->rLambFront;
@@ -565,7 +572,7 @@ SDevalBSDF(SDValue *sv, const FVECT outVec, const FVEC
 	i = (sdf != NULL) ? sdf->ncomp : 0;
 	while (i-- > 0) {
 		nch = (*sdf->comp[i].func->getBSDFs)(coef, outVec, inVec,
-							sdf->comp[i].dist);
+							&sdf->comp[i]);
 		while (nch-- > 0) {
 			c_cmix(&sv->spec, sv->cieY, &sv->spec,
 					coef[nch], &sdf->comp[i].cspec[nch]);
@@ -589,7 +596,7 @@ SDdirectHemi(const FVECT inVec, int sflags, const SDDa
 	if ((inVec == NULL) | (sd == NULL))
 		return .0;
 					/* gather diffuse components */
-	if (inVec[2] > .0) {
+	if (inVec[2] > 0) {
 		hsum = sd->rLambFront.cieY;
 		rdf = sd->rf;
 	} else /* !inFront */ {
@@ -620,10 +627,10 @@ SDdirectHemi(const FVECT inVec, int sflags, const SDDa
 
 /* Sample BSDF direction based on the given random variable */
 SDError
-SDsampBSDF(SDValue *sv, FVECT outVec, const FVECT inVec,
-			double randX, int sflags, const SDData *sd)
+SDsampBSDF(SDValue *sv, FVECT ioVec, double randX, int sflags, const SDData *sd)
 {
 	SDError		ec;
+	FVECT		inVec;
 	int		inFront;
 	SDSpectralDF	*rdf;
 	double		rdiff;
@@ -632,11 +639,12 @@ SDsampBSDF(SDValue *sv, FVECT outVec, const FVECT inVe
 	SDComponent	*sdc;
 	const SDCDst	**cdarr = NULL;
 					/* check arguments */
-	if ((sv == NULL) | (outVec == NULL) | (inVec == NULL) | (sd == NULL) |
-			(randX < .0) | (randX >= 1.))
+	if ((sv == NULL) | (ioVec == NULL) | (sd == NULL) |
+			(randX < 0) | (randX >= 1.))
 		return SDEargument;
 					/* whose side are we on? */
-	inFront = (inVec[2] > .0);
+	VCOPY(inVec, ioVec);
+	inFront = (inVec[2] > 0);
 					/* remember diffuse portions */
 	if (inFront) {
 		*sv = sd->rLambFront;
@@ -676,14 +684,14 @@ SDsampBSDF(SDValue *sv, FVECT outVec, const FVECT inVe
 	}
 	if (sv->cieY <= 1e-7) {		/* anything to sample? */
 		sv->cieY = .0;
-		memset(outVec, 0, 3*sizeof(double));
+		memset(ioVec, 0, 3*sizeof(double));
 		return SDEnone;
 	}
 					/* scale random variable */
 	randX *= sv->cieY;
 					/* diffuse reflection? */
 	if (randX < rdiff) {
-		SDdiffuseSamp(outVec, inFront, randX/rdiff);
+		SDdiffuseSamp(ioVec, inFront, randX/rdiff);
 		goto done;
 	}
 	randX -= rdiff;
@@ -691,7 +699,7 @@ SDsampBSDF(SDValue *sv, FVECT outVec, const FVECT inVe
 	if ((sflags & SDsampDf+SDsampT) == SDsampDf+SDsampT) {
 		if (randX < sd->tLamb.cieY) {
 			sv->spec = sd->tLamb.spec;
-			SDdiffuseSamp(outVec, !inFront, randX/sd->tLamb.cieY);
+			SDdiffuseSamp(ioVec, !inFront, randX/sd->tLamb.cieY);
 			goto done;
 		}
 		randX -= sd->tLamb.cieY;
@@ -703,11 +711,11 @@ SDsampBSDF(SDValue *sv, FVECT outVec, const FVECT inVe
 		return SDEinternal;
 					/* compute sample direction */
 	sdc = (i < nr) ? &rdf->comp[i] : &sd->tf->comp[i-nr];
-	ec = (*sdc->func->sampCDist)(outVec, randX/cdarr[i]->cTotal, cdarr[i]);
+	ec = (*sdc->func->sampCDist)(ioVec, randX/cdarr[i]->cTotal, cdarr[i]);
 	if (ec)
 		return ec;
 					/* compute color */
-	j = (*sdc->func->getBSDFs)(coef, outVec, inVec, sdc->dist);
+	j = (*sdc->func->getBSDFs)(coef, ioVec, inVec, sdc);
 	if (j <= 0) {
 		sprintf(SDerrorDetail, "BSDF \"%s\" sampling value error",
 				sd->name);
@@ -734,10 +742,10 @@ SDcompXform(RREAL vMtx[3][3], const FVECT sNrm, const 
 	if ((vMtx == NULL) | (sNrm == NULL) | (uVec == NULL))
 		return SDEargument;
 	VCOPY(vMtx[2], sNrm);
-	if (normalize(vMtx[2]) == .0)
+	if (normalize(vMtx[2]) == 0)
 		return SDEargument;
 	fcross(vMtx[0], uVec, vMtx[2]);
-	if (normalize(vMtx[0]) == .0)
+	if (normalize(vMtx[0]) == 0)
 		return SDEargument;
 	fcross(vMtx[1], vMtx[2], vMtx[0]);
 	return SDEnone;
@@ -757,7 +765,7 @@ SDinvXform(RREAL iMtx[3][3], RREAL vMtx[3][3])
 	mTmp[0][1] = vMtx[2][1]*vMtx[0][2] - vMtx[2][2]*vMtx[0][1];
 	mTmp[0][2] = vMtx[1][2]*vMtx[0][1] - vMtx[1][1]*vMtx[0][2];
 	d = vMtx[0][0]*mTmp[0][0] + vMtx[1][0]*mTmp[0][1] + vMtx[2][0]*mTmp[0][2];
-	if (d == .0) {
+	if (d == 0) {
 		strcpy(SDerrorDetail, "Zero determinant in matrix inversion");
 		return SDEargument;
 	}
@@ -784,12 +792,12 @@ SDmapDir(FVECT resVec, RREAL vMtx[3][3], const FVECT i
 	if (vMtx == NULL) {		/* assume they just want to normalize */
 		if (resVec != inpVec)
 			VCOPY(resVec, inpVec);
-		return (normalize(resVec) > .0) ? SDEnone : SDEargument;
+		return (normalize(resVec) > 0) ? SDEnone : SDEargument;
 	}
 	vTmp[0] = DOT(vMtx[0], inpVec);
 	vTmp[1] = DOT(vMtx[1], inpVec);
 	vTmp[2] = DOT(vMtx[2], inpVec);
-	if (normalize(vTmp) == .0)
+	if (normalize(vTmp) == 0)
 		return SDEargument;
 	VCOPY(resVec, vTmp);
 	return SDEnone;
@@ -861,7 +869,7 @@ static int	nabases = 3;	/* current number of defined b
 static int
 fequal(double a, double b)
 {
-	if (b != .0)
+	if (b != 0)
 		a = a/b - 1.;
 	return((a <= 1e-6) & (a >= -1e-6));
 }
@@ -1179,7 +1187,7 @@ check_bsdf_data(	/* check that BSDF data is sane */
 	hemi_total = .0;
 	for (i = dp->ninc; i--; ) {
 		dom = getBSDF_incohm(dp,i);
-		if (dom <= .0) {
+		if (dom <= 0) {
 			error(WARNING, "zero/negative incoming solid angle");
 			continue;
 		}
@@ -1202,7 +1210,7 @@ check_bsdf_data(	/* check that BSDF data is sane */
 	hemi_total = .0;
 	for (o = dp->nout; o--; ) {
 		dom = getBSDF_outohm(dp,o);
-		if (dom <= .0) {
+		if (dom <= 0) {
 			error(WARNING, "zero/negative outgoing solid angle");
 			continue;
 		}
@@ -1226,7 +1234,7 @@ check_bsdf_data(	/* check that BSDF data is sane */
 		hemi_total = .0;
 		for (o = dp->nout; o--; ) {
 			double	f = BSDF_value(dp,i,o);
-			if (f >= .0)
+			if (f >= 0)
 				hemi_total += f*omega_oarr[o];
 			else {
 				nneg += (f < -FTINY);