--- ray/src/util/dctimestep.c	2013/01/11 05:07:47	2.24
+++ ray/src/util/dctimestep.c	2013/01/19 22:30:55	2.27
@@ -1,5 +1,5 @@
 #ifndef lint
-static const char RCSid[] = "$Id: dctimestep.c,v 2.24 2013/01/11 05:07:47 greg Exp $";
+static const char RCSid[] = "$Id: dctimestep.c,v 2.27 2013/01/19 22:30:55 greg Exp $";
 #endif
 /*
  * Compute time-step result using Daylight Coefficient method.
@@ -107,8 +107,8 @@ getDTfromHeader(FILE *fp)
 static CMATRIX *
 cm_load(const char *fname, int nrows, int ncols, int dtype)
 {
-	CMATRIX	*cm;
 	FILE	*fp = stdin;
+	CMATRIX	*cm;
 
 	if (ncols <= 0)
 		error(USER, "Non-positive number of columns");
@@ -118,8 +118,11 @@ cm_load(const char *fname, int nrows, int ncols, int d
 		sprintf(errmsg, "cannot open file '%s'", fname);
 		error(SYSTEM, errmsg);
 	}
+#ifdef getc_unlocked
+	flockfile(fp);
+#endif
 	if (dtype != DTascii)
-		SET_FILE_BINARY(fp);
+		SET_FILE_BINARY(fp);		/* doesn't really work */
 	if (dtype == DTfromHeader)
 		dtype = getDTfromHeader(fp);
 	switch (dtype) {
@@ -158,7 +161,7 @@ cm_load(const char *fname, int nrows, int ncols, int d
 		cm = cm_alloc(guessrows, ncols);
 	} else
 		cm = cm_alloc(nrows, ncols);
-	if (cm == NULL)
+	if (cm == NULL)					/* XXX never happens */
 		return(NULL);
 	if (dtype == DTascii) {				/* read text file */
 		int	maxrow = (nrows > 0 ? nrows : 32000);
@@ -233,7 +236,7 @@ cm_load(const char *fname, int nrows, int ncols, int d
 				cvp += 3;
 			}
 		}
-		if (getc(fp) != EOF) {
+		if (fgetc(fp) != EOF) {
 				sprintf(errmsg,
 				"unexpected data at end of binary file %s",
 						fname);
@@ -242,10 +245,13 @@ cm_load(const char *fname, int nrows, int ncols, int d
 	}
 	if (fp != stdin)
 		fclose(fp);
+#ifdef getc_unlocked
+	else
+		funlockfile(fp);
+#endif
 	return(cm);
 EOFerror:
-	sprintf(errmsg, "unexpected EOF reading %s",
-			fname);
+	sprintf(errmsg, "unexpected EOF reading %s", fname);
 	error(USER, errmsg);
 not_handled:
 	error(INTERNAL, "unhandled data size or length in cm_load()");
@@ -260,7 +266,7 @@ cm_column(const CMATRIX *cm, int c)
 	int	dr;
 
 	if ((c < 0) | (c >= cm->ncols))
-		return(NULL);
+		error(INTERNAL, "column requested outside matrix");
 	cvr = cm_alloc(cm->nrows, 1);
 	if (cvr == NULL)
 		return(NULL);
@@ -299,6 +305,7 @@ cm_scale(const CMATRIX *cm1, const COLOR sca)
 static CMATRIX *
 cm_multiply(const CMATRIX *cm1, const CMATRIX *cm2)
 {
+	char	*rowcheck=NULL, *colcheck=NULL;
 	CMATRIX	*cmr;
 	int	dr, dc, i;
 
@@ -307,10 +314,32 @@ cm_multiply(const CMATRIX *cm1, const CMATRIX *cm2)
 	cmr = cm_alloc(cm1->nrows, cm2->ncols);
 	if (cmr == NULL)
 		return(NULL);
+				/* optimization: check for zero rows & cols */
+	if (((cm1->nrows > 5) | (cm2->ncols > 5)) & (cm1->ncols > 5)) {
+		static const COLOR	czero;
+		rowcheck = (char *)calloc(cmr->nrows, 1);
+		for (dr = cm1->nrows*(rowcheck != NULL); dr--; )
+		    for (dc = cm1->ncols; dc--; )
+			if (memcmp(cm_lval(cm1,dr,dc), czero, sizeof(COLOR))) {
+				rowcheck[dr] = 1;
+				break;
+			}
+		colcheck = (char *)calloc(cmr->ncols, 1);
+		for (dc = cm2->ncols*(colcheck != NULL); dc--; )
+		    for (dr = cm2->nrows; dr--; )
+			if (memcmp(cm_lval(cm2,dr,dc), czero, sizeof(COLOR))) {
+				colcheck[dc] = 1;
+				break;
+			}
+	}
 	for (dr = 0; dr < cmr->nrows; dr++)
 	    for (dc = 0; dc < cmr->ncols; dc++) {
 		COLORV	*dp = cm_lval(cmr,dr,dc);
 		dp[0] = dp[1] = dp[2] = 0;
+		if (rowcheck != NULL && !rowcheck[dr])
+			continue;
+		if (colcheck != NULL && !colcheck[dc])
+			continue;
 		for (i = 0; i < cm1->ncols; i++) {
 		    const COLORV	*cp1 = cm_lval(cm1,dr,i);
 		    const COLORV	*cp2 = cm_lval(cm2,i,dc);
@@ -319,6 +348,8 @@ cm_multiply(const CMATRIX *cm1, const CMATRIX *cm2)
 		    dp[2] += cp1[2] * cp2[2];
 		}
 	    }
+	if (rowcheck != NULL) free(rowcheck);
+	if (colcheck != NULL) free(colcheck);
 	return(cmr);
 }
 
@@ -587,6 +618,7 @@ hasNumberFormat(const char *s)
 int
 main(int argc, char *argv[])
 {
+	int		skyfmt = DTascii;
 	int		nsteps = 1;
 	char		*ofspec = NULL;
 	FILE		*ofp = stdout;
@@ -596,8 +628,8 @@ main(int argc, char *argv[])
 
 	progname = argv[0];
 					/* get options */
-	for (a = 1; a < argc-1 && argv[a][0] == '-'; a++)
-		switch (argv[0][1]) {
+	for (a = 1; a < argc && argv[a][0] == '-'; a++)
+		switch (argv[a][1]) {
 		case 'n':
 			nsteps = atoi(argv[++a]);
 			if (nsteps <= 0)
@@ -606,6 +638,21 @@ main(int argc, char *argv[])
 		case 'o':
 			ofspec = argv[++a];
 			break;
+		case 'i':
+			switch (argv[a][2]) {
+			case 'f':
+				skyfmt = DTfloat;
+				break;
+			case 'd':
+				skyfmt = DTdouble;
+				break;
+			case 'a':
+				skyfmt = DTascii;
+				break;
+			default:
+				goto userr;
+			}
+			break;
 		default:
 			goto userr;
 		}
@@ -616,7 +663,7 @@ main(int argc, char *argv[])
 		CMATRIX	*smtx, *Dmat, *Tmat, *imtx;
 		COLOR	tLamb;
 						/* get sky vector/matrix */
-		smtx = cm_load(argv[a+3], 0, nsteps, DTascii);
+		smtx = cm_load(argv[a+3], 0, nsteps, skyfmt);
 						/* load BSDF */
 		Tmat = cm_loadBSDF(argv[a+1], tLamb);
 						/* load Daylight matrix */
@@ -633,7 +680,7 @@ main(int argc, char *argv[])
 		}
 		cm_free(imtx);
 	} else {				/* sky vector/matrix only */
-		cmtx = cm_load(argv[a+1], 0, nsteps, DTascii);
+		cmtx = cm_load(argv[a+1], 0, nsteps, skyfmt);
 	}
 						/* prepare output stream */
 	if ((ofspec != NULL) & (nsteps == 1) && hasNumberFormat(ofspec)) {
@@ -721,9 +768,9 @@ main(int argc, char *argv[])
 	cm_free(cmtx);
 	return(0);
 userr:
-	fprintf(stderr, "Usage: %s [-n nsteps][-o ospec] DCspec [skyf]\n",
+	fprintf(stderr, "Usage: %s [-n nsteps][-o ospec][-i{f|d}] DCspec [skyf]\n",
 				progname);
-	fprintf(stderr, "   or: %s [-n nsteps][-o ospec] Vspec Tbsdf.xml Dmat.dat [skyf]\n",
+	fprintf(stderr, "   or: %s [-n nsteps][-o ospec][-i{f|d}] Vspec Tbsdf.xml Dmat.dat [skyf]\n",
 				progname);
 	return(1);
 }