--- ray/src/util/dctimestep.c	2011/04/08 18:13:48	2.15
+++ ray/src/util/dctimestep.c	2012/10/16 15:53:04	2.21
@@ -1,5 +1,5 @@
 #ifndef lint
-static const char RCSid[] = "$Id: dctimestep.c,v 2.15 2011/04/08 18:13:48 greg Exp $";
+static const char RCSid[] = "$Id: dctimestep.c,v 2.21 2012/10/16 15:53:04 greg Exp $";
 #endif
 /*
  * Compute time-step result using Daylight Coefficient method.
@@ -252,8 +252,29 @@ not_handled:
 	return(NULL);	/* gratis return */
 }
 
-/* Multiply two matrices (or a matrix and a vector) and allocate the result*/
+/* Scale a matrix by a single value */
 static CMATRIX *
+cm_scale(const CMATRIX *cm1, const COLOR sca)
+{
+	CMATRIX	*cmr;
+	int	dr, dc;
+
+	cmr = cm_alloc(cm1->nrows, cm1->ncols);
+	if (cmr == NULL)
+		return(NULL);
+	for (dr = 0; dr < cmr->nrows; dr++)
+	    for (dc = 0; dc < cmr->ncols; dc++) {
+	        const COLORV	*sp = cm_lval(cm1,dr,dc);
+		COLORV		*dp = cm_lval(cmr,dr,dc);
+		dp[0] = sp[0] * sca[0];
+		dp[1] = sp[1] * sca[1];
+		dp[2] = sp[2] * sca[2];
+	    }
+	return(cmr);
+}
+
+/* Multiply two matrices (or a matrix and a vector) and allocate the result */
+static CMATRIX *
 cm_multiply(const CMATRIX *cm1, const CMATRIX *cm2)
 {
 	CMATRIX	*cmr;
@@ -297,66 +318,29 @@ cm_print(const CMATRIX *cm, FILE *fp)
 static CMATRIX *
 cm_bsdf(const COLOR bsdfLamb, const COLOR specCol, const SDMat *bsdf)
 {
-	CMATRIX	*cm;
+	CMATRIX	*cm = cm_alloc(bsdf->nout, bsdf->ninc);
 	int	nbadohm = 0;
 	int	nneg = 0;
-	FVECT	v;
-	float	dom;
-	int	doforward;
 	int	r, c;
+					/* reciprocity is "transparent" */
+	for (c = 0; c < cm->ncols; c++) {
+		const double	dom = mBSDF_incohm(bsdf,c);
+					/* projected solid angle */
+		nbadohm += (dom <= 0);
 
-	mBSDF_incvec(v, bsdf, .5);	/* check BTDF orientation */
-
-	if (v[2] < 0) {			/* incident from outside */
-		cm = cm_alloc(bsdf->nout, bsdf->ninc);
-		for (c = 0; c < cm->ncols; c++) {
-			dom = mBSDF_incohm(bsdf,c);
-			
-			nbadohm += (dom <= 0);
-
-			if (!mBSDF_incvec(v,bsdf,c+.5) || v[2] > 0)
-				error(USER, "illegal incoming BTDF direction");
-			dom *= -v[2];
-
-			for (r = 0; r < cm->nrows; r++) {
-				float	f = mBSDF_value(bsdf,c,r);
-				COLORV	*mp = cm_lval(cm,r,c);
-
-				if ((f <= 0) | (dom <= 0)) {
-					nneg += (f < -FTINY);
-					f = .0f;
-				}
-				copycolor(mp, specCol);
-				f *= dom;
-				scalecolor(mp, f);
-				addcolor(mp, bsdfLamb);
+		for (r = 0; r < cm->nrows; r++) {
+			float	f = mBSDF_value(bsdf,c,r);
+			COLORV	*mp = cm_lval(cm,r,c);
+					/* check BSDF value */
+			if ((f <= 0) | (dom <= 0)) {
+				nneg += (f < -FTINY);
+				f = .0f;
 			}
+			copycolor(mp, specCol);
+			scalecolor(mp, f);
+			addcolor(mp, bsdfLamb);
+			scalecolor(mp, dom);
 		}
-	} else {			/* rely on reciprocity */
-		cm = cm_alloc(bsdf->ninc, bsdf->nout);
-		for (c = 0; c < cm->ncols; c++) {
-			dom = mBSDF_outohm(bsdf,c);
-			
-			nbadohm += (dom <= 0);
-
-			if (!mBSDF_outvec(v,bsdf,c+.5) || v[2] > 0)
-				error(USER, "illegal outgoing BTDF direction");
-			dom *= -v[2];
-
-			for (r = 0; r < cm->nrows; r++) {
-				float	f = mBSDF_value(bsdf,r,c);
-				COLORV	*mp = cm_lval(cm,r,c);
-
-				if ((f <= 0) | (dom <= 0)) {
-					nneg += (f < -FTINY);
-					f = .0f;
-				}
-				copycolor(mp, specCol);
-				f *= dom;
-				scalecolor(mp, f);
-				addcolor(mp, bsdfLamb);
-			}
-		}
 	}
 	if (nneg | nbadohm) {
 		sprintf(errmsg,
@@ -369,33 +353,41 @@ cm_bsdf(const COLOR bsdfLamb, const COLOR specCol, con
 
 /* Load and convert a matrix BSDF from the given XML file */
 static CMATRIX *
-cm_loadBSDF(char *fname)
+cm_loadBSDF(char *fname, COLOR cLamb)
 {
-	CMATRIX	*Tmat;
-	char	*fpath;
-	SDError	ec;
-	SDData	myBSDF;
-	COLOR	bsdfLamb, specCol;
+	CMATRIX		*Tmat;
+	char		*fpath;
+	SDError		ec;
+	SDData		myBSDF;
+	SDSpectralDF	*tdf;
+	COLOR		bsdfLamb, specCol;
 					/* find path to BSDF file */
 	fpath = getpath(fname, getrlibpath(), R_OK);
 	if (fpath == NULL) {
 		sprintf(errmsg, "cannot find BSDF file '%s'", fname);
 		error(USER, errmsg);
 	}
-	SDclipName(myBSDF.name, fname);
-					/* load XML and check type */
+	SDclearBSDF(&myBSDF, fname);	/* load XML and check type */
 	ec = SDloadFile(&myBSDF, fpath);
 	if (ec)
 		error(USER, transSDError(ec));
-	if (myBSDF.tf == NULL || myBSDF.tf->ncomp != 1 ||
-			myBSDF.tf->comp[0].func != &SDhandleMtx) {
+	ccy2rgb(&myBSDF.tLamb.spec, myBSDF.tLamb.cieY/PI, bsdfLamb);
+	tdf = (myBSDF.tf != NULL) ? myBSDF.tf : myBSDF.tb;
+	if (tdf == NULL) {		/* no non-Lambertian transmission? */
+		if (cLamb != NULL)
+			copycolor(cLamb, bsdfLamb);
+		SDfreeBSDF(&myBSDF);
+		return(NULL);
+	}
+	if (tdf->ncomp != 1 || tdf->comp[0].func != &SDhandleMtx) {
 		sprintf(errmsg, "unsupported BSDF '%s'", fpath);
 		error(USER, errmsg);
 	}
 					/* convert BTDF to matrix */
-	ccy2rgb(&myBSDF.tLamb.spec, myBSDF.tLamb.cieY/PI, bsdfLamb);
-	ccy2rgb(&myBSDF.tf->comp[0].cspec[0], 1., specCol);
-	Tmat = cm_bsdf(bsdfLamb, specCol, (SDMat *)myBSDF.tf->comp[0].dist);
+	ccy2rgb(&tdf->comp[0].cspec[0], 1., specCol);
+	Tmat = cm_bsdf(bsdfLamb, specCol, (SDMat *)tdf->comp[0].dist);
+	if (cLamb != NULL)		/* Lambertian is included */
+		setcolor(cLamb, .0, .0, .0);
 					/* free BSDF and return */
 	SDfreeBSDF(&myBSDF);
 	return(Tmat);
@@ -420,7 +412,8 @@ sum_images(const char *fspec, const CMATRIX *cv, FILE 
 		int		dt, xr, yr;
 		COLORV		*psp;
 							/* check for zero */
-		if ((scv[RED] == 0) & (scv[GRN] == 0) & (scv[BLU] == 0))
+		if ((scv[RED] == 0) & (scv[GRN] == 0) & (scv[BLU] == 0) &&
+				(myDT != DTfromHeader) | (i < cv->nrows-1))
 			continue;
 							/* open next picture */
 		sprintf(fname, fspec, i);
@@ -482,22 +475,28 @@ sum_images(const char *fspec, const CMATRIX *cv, FILE 
 static int
 hasNumberFormat(const char *s)
 {
-	while (*s && *s != '%')
-		s++;
-	if (!*s)
-		return(0);
-	do
-		++s;
-	while (isdigit(*s));
-
-	return((*s == 'd') | (*s == 'i') | (*s == 'o') |
-			(*s == 'x') | (*s == 'X'));
+	while (*s) {
+		while (*s != '%')
+			if (!*s++)
+				return(0);
+		if (*++s == '%') {		/* ignore "%%" */
+			++s;
+			continue;
+		}
+		while (isdigit(*s))		/* field length */
+			++s;
+						/* field we'll use? */
+		if ((*s == 'd') | (*s == 'i') | (*s == 'o') |
+					(*s == 'x') | (*s == 'X'))
+			return(1);
+	}
+	return(0);				/* didn't find one */
 }
 
 int
 main(int argc, char *argv[])
 {
-	CMATRIX			*cvec;
+	CMATRIX			*cvec;		/* component vector result */
 
 	progname = argv[0];
 
@@ -510,17 +509,24 @@ main(int argc, char *argv[])
 
 	if (argc > 3) {				/* VTDs expression */
 		CMATRIX	*svec, *Dmat, *Tmat, *ivec;
+		COLOR	tLamb;
 						/* get sky vector */
 		svec = cm_load(argv[4], 0, 1, DTascii);
 						/* load BSDF */
-		Tmat = cm_loadBSDF(argv[2]);
+		Tmat = cm_loadBSDF(argv[2], tLamb);
 						/* load Daylight matrix */
-		Dmat = cm_load(argv[3], Tmat->ncols, svec->nrows, DTfromHeader);
+		Dmat = cm_load(argv[3], Tmat==NULL ? 0 : Tmat->ncols,
+					svec->nrows, DTfromHeader);
 						/* multiply vector through */
 		ivec = cm_multiply(Dmat, svec);
 		cm_free(Dmat); cm_free(svec);
-		cvec = cm_multiply(Tmat, ivec);	/* cvec = component vector */
-		cm_free(Tmat); cm_free(ivec);
+		if (Tmat == NULL) {		/* diffuse only */
+			cvec = cm_scale(ivec, tLamb);
+		} else {			/* else apply BTDF matrix */
+			cvec = cm_multiply(Tmat, ivec);
+			cm_free(Tmat); 
+		}
+		cm_free(ivec);
 	} else {				/* else just use sky vector */
 		cvec = cm_load(argv[2], 0, 1, DTascii);
 	}