--- ray/src/cv/bsdf2rad.c	2017/04/12 03:35:01	2.26
+++ ray/src/cv/bsdf2rad.c	2021/12/07 23:55:02	2.39
@@ -1,21 +1,22 @@
 #ifndef lint
-static const char RCSid[] = "$Id: bsdf2rad.c,v 2.26 2017/04/12 03:35:01 greg Exp $";
+static const char RCSid[] = "$Id: bsdf2rad.c,v 2.39 2021/12/07 23:55:02 greg Exp $";
 #endif
 /*
  *  Plot 3-D BSDF output based on scattering interpolant or XML representation
  */
 
-#include <stdio.h>
-#include <string.h>
 #include <stdlib.h>
+#include "rtio.h"
 #include "paths.h"
 #include "rtmath.h"
-#include "resolu.h"
 #include "bsdfrep.h"
 
+#ifndef NINCIDENT
 #define NINCIDENT	37		/* number of samples/hemisphere */
-
+#endif
+#ifndef GRIDSTEP
 #define	GRIDSTEP	2		/* our grid step size */
+#endif
 #define SAMPRES		(GRIDRES/GRIDSTEP)
 
 int	front_comp = 0;			/* front component flags (SDsamp*) */
@@ -24,7 +25,7 @@ double	overall_min = 1./PI;		/* overall minimum BSDF v
 double	min_log10;			/* smallest log10 value for plotting */
 double	overall_max = .0;		/* overall maximum BSDF value */
 
-char	ourTempDir[TEMPLEN] = "";	/* our temporary directory */
+char	ourTempDir[TEMPLEN+1] = "";	/* our temporary directory */
 
 const char	frpref[] = "rf";
 const char	ftpref[] = "tf";
@@ -40,10 +41,9 @@ const double	sph_xoffset = 15.;
 #define bsdf_rad	(sph_rad*.25)
 #define arrow_rad	(bsdf_rad*.015)
 
-#define	FEQ(a,b)	((a)-(b) <= 1e-7 && (b)-(a) <= 1e-7)
+#define	set_minlog()	overall_min = (overall_min < 1e-5) ? 1e-5 : overall_min; \
+				min_log10 = log10(overall_min) - .1
 
-#define	set_minlog()	(min_log10 = log10(overall_min + 1e-5) - .1)
-
 char	*progname;
 
 /* Get Fibonacci sphere vector (0 to NINCIDENT-1) */
@@ -76,7 +76,7 @@ cvt_sposition(FVECT sp, const FVECT iv, int inc_side)
 static char *
 tfile_name(const char *prefix, const char *suffix, int i)
 {
-	static char	buf[128];
+	static char	buf[256];
 
 	if (!ourTempDir[0]) {		/* create temporary directory */
 		mktemp(strcpy(ourTempDir,TEMPLATE));
@@ -144,8 +144,8 @@ plotBSDF(const char *fname, const FVECT ivec, int dfl,
 		SDValue	sval;
 		double	bsdf;
 		ovec_from_pos(ovec, i*GRIDSTEP, j*GRIDSTEP);
-		if (SDreportError(SDevalBSDF(&sval, ovec,
-						ivec, sd), stderr))
+		if (SDreportError(SDevalBSDF(&sval, ivec,
+						ovec, sd), stderr))
 			return(0);
 		if (sval.cieY > overall_max)
 			overall_max = sval.cieY;
@@ -274,6 +274,7 @@ put_mirror_arrow(const FVECT origin, const FVECT nrm)
 {
 	const double	arrow_len = 1.2*bsdf_rad;
 	const double	tip_len = 0.2*bsdf_rad;
+	static int	cnt = 1;
 	FVECT		refl;
 	int		i;
 
@@ -281,20 +282,20 @@ put_mirror_arrow(const FVECT origin, const FVECT nrm)
 	refl[1] = 2.*nrm[2]*nrm[1];
 	refl[2] = 2.*nrm[2]*nrm[2] - 1.;
 
-	printf("\n# Mirror arrow\n");
-	printf("\narrow_mat cylinder inc_dir\n0\n0\n7");
+	printf("\n# Mirror arrow #%d\n", cnt);
+	printf("\nshaft_mat cylinder inc_dir%d\n0\n0\n7", cnt);
 	printf("\n\t%f %f %f\n\t%f %f %f\n\t%f\n",
 			origin[0], origin[1], origin[2]+arrow_len,
 			origin[0], origin[1], origin[2],
 			arrow_rad);
-	printf("\narrow_mat cylinder mir_dir\n0\n0\n7");
+	printf("\nshaft_mat cylinder mir_dir%d\n0\n0\n7", cnt);
 	printf("\n\t%f %f %f\n\t%f %f %f\n\t%f\n",
 			origin[0], origin[1], origin[2],
 			origin[0] + arrow_len*refl[0],
 			origin[1] + arrow_len*refl[1],
 			origin[2] + arrow_len*refl[2],
 			arrow_rad);
-	printf("\narrow_mat cone mir_tip\n0\n0\n8");
+	printf("\ntip_mat cone mir_tip%d\n0\n0\n8", cnt);
 	printf("\n\t%f %f %f\n\t%f %f %f\n\t%f 0\n",
 			origin[0] + (arrow_len-.5*tip_len)*refl[0],
 			origin[1] + (arrow_len-.5*tip_len)*refl[1],
@@ -303,6 +304,7 @@ put_mirror_arrow(const FVECT origin, const FVECT nrm)
 			origin[1] + (arrow_len+.5*tip_len)*refl[1],
 			origin[2] + (arrow_len+.5*tip_len)*refl[2],
 			2.*arrow_rad);
+	++cnt;
 }
 
 /* Put out transmitted direction arrow for the given incident vector */
@@ -311,19 +313,21 @@ put_trans_arrow(const FVECT origin)
 {
 	const double	arrow_len = 1.2*bsdf_rad;
 	const double	tip_len = 0.2*bsdf_rad;
+	static int	cnt = 1;
 	int		i;
 
-	printf("\n# Transmission arrow\n");
-	printf("\narrow_mat cylinder trans_dir\n0\n0\n7");
+	printf("\n# Transmission arrow #%d\n", cnt);
+	printf("\nshaft_mat cylinder trans_dir%d\n0\n0\n7", cnt);
 	printf("\n\t%f %f %f\n\t%f %f %f\n\t%f\n",
 			origin[0], origin[1], origin[2],
 			origin[0], origin[1], origin[2]-arrow_len,
 			arrow_rad);
-	printf("\narrow_mat cone trans_tip\n0\n0\n8");
+	printf("\ntip_mat cone trans_tip%d\n0\n0\n8", cnt);
 	printf("\n\t%f %f %f\n\t%f %f %f\n\t%f 0\n",
 			origin[0], origin[1], origin[2]-arrow_len+.5*tip_len,
 			origin[0], origin[1], origin[2]-arrow_len-.5*tip_len,
-			2.*arrow_rad);	
+			2.*arrow_rad);
+	++cnt;
 }
 
 /* Compute rotation (x,y,z) => (xp,yp,zp) */
@@ -342,19 +346,19 @@ addrot(char *xf, const FVECT xp, const FVECT yp, const
 		return(4);
 	}
 	theta = atan2(yp[2], zp[2]);
-	if (!FEQ(theta,0.0)) {
+	if (!FABSEQ(theta,0.0)) {
 		sprintf(xf, " -rx %f", theta*(180./PI));
 		while (*xf) ++xf;
 		n += 2;
 	}
 	theta = Asin(-xp[2]);
-	if (!FEQ(theta,0.0)) {
+	if (!FABSEQ(theta,0.0)) {
 		sprintf(xf, " -ry %f", theta*(180./PI));
 		while (*xf) ++xf;
 		n += 2;
 	}
 	theta = atan2(xp[1], xp[0]);
-	if (!FEQ(theta,0.0)) {
+	if (!FABSEQ(theta,0.0)) {
 		sprintf(xf, " -rz %f", theta*(180./PI));
 		/* while (*xf) ++xf; */
 		n += 2;
@@ -378,8 +382,8 @@ put_BSDFs(void)
 	printf("\n# Gensurf output corresponding to %d incident directions\n",
 			NINCIDENT);
 
-	printf("\nvoid glow arrow_glow\n0\n0\n4 1 0 1 0\n");
-	printf("\nvoid mixfunc arrow_mat\n4 arrow_glow void 0.25 .\n0\n0\n");
+	printf("\nvoid glow tip_mat\n0\n0\n4 1 0 1 0\n");
+	printf("\nvoid mixfunc shaft_mat\n4 tip_mat void 0.25 .\n0\n0\n");
 
 	for (i = 0; i < NINCIDENT; i++) {
 		get_ivector(ivec, i);
@@ -470,6 +474,7 @@ put_matBSDF(const char *XMLfile)
 	}
 	switch (XMLfile[0]) {		/* avoid RAYPATH search */
 	case '.':
+	case '~':
 	CASEDIRSEP:
 		curdir = "";
 		break;
@@ -511,7 +516,7 @@ put_hemispheres(void)
 	if (front_comp) {
 		printf(
 "\n!genrev %s Front \"R*sin(A*t)\" \"R*cos(A*t)\" %d -e \"R:%g;A:%f\" -s | xform -t %g 0 0\n",
-				sph_fmat, nsegs, sph_rad, 0.495*PI, sph_xoffset);
+				sph_fmat, nsegs, sph_rad, 0.5*PI, sph_xoffset);
 		printf("\nvoid brighttext front_text\n3 helvet.fnt . FRONT\n0\n");
 		printf("12\n\t%f %f 0\n\t%f 0 0\n\t0 %f 0\n\t.01 1 -.1\n",
 				-.22*sph_rad + sph_xoffset, -1.4*sph_rad,
@@ -527,7 +532,7 @@ put_hemispheres(void)
 	if (back_comp) {
 		printf(
 "\n!genrev %s Back \"R*cos(A*t)\" \"R*sin(A*t)\" %d -e \"R:%g;A:%f\" -s | xform -t %g 0 0\n",
-				sph_bmat, nsegs, sph_rad, 0.495*PI, -sph_xoffset);
+				sph_bmat, nsegs, sph_rad, 0.5*PI, -sph_xoffset);
 		printf("\nvoid brighttext back_text\n3 helvet.fnt . BACK\n0\n");
 		printf("12\n\t%f %f 0\n\t%f 0 0\n\t0 %f 0\n\t.01 1 -.1\n",
 				-.22*sph_rad - sph_xoffset, -1.4*sph_rad,
@@ -667,7 +672,7 @@ convert_mgf(const char *mgfdata)
 static int
 rbf_headline(char *s, void *p)
 {
-	char	fmt[64];
+	char	fmt[MAXFMTLEN];
 
 	if (formatval(fmt, s)) {
 		if (strcmp(fmt, BSDFREP_FMT))
@@ -701,27 +706,35 @@ int
 main(int argc, char *argv[])
 {
 	int	inpXML = -1;
+	double	myLim[2];
 	SDData	myBSDF;
-	int	n;
+	int	a, n;
 						/* check arguments */
 	progname = argv[0];
-	if (argc > 1 && (n = strlen(argv[1])-4) > 0) {
-		if (!strcasecmp(argv[1]+n, ".xml"))
+	a = 1;
+	myLim[0] = -1; myLim[1] = -2;		/* specified BSDF range? */
+	if (argc > a+3 && argv[a][0] == '-' && argv[a][1] == 'r') {
+		myLim[0] = atof(argv[++a]);
+		myLim[1] = atof(argv[++a]);
+		++a;
+	}
+	if (argc > a && (n = strlen(argv[a])-4) > 0) {
+		if (!strcasecmp(argv[a]+n, ".xml"))
 			inpXML = 1;
-		else if (!strcasecmp(argv[1]+n, ".sir"))
+		else if (!strcasecmp(argv[a]+n, ".sir"))
 			inpXML = 0;
 	}
-	if (inpXML < 0 || inpXML & (argc > 2)) {
-		fprintf(stderr, "Usage: %s bsdf.xml > output.rad\n", progname);
-		fprintf(stderr, "   Or: %s hemi1.sir hemi2.sir .. > output.rad\n", progname);
+	if (inpXML < 0 || inpXML & (argc > a+1)) {
+		fprintf(stderr, "Usage: %s [-r min max] bsdf.xml > output.rad\n", progname);
+		fprintf(stderr, "   Or: %s [-r min max] hemi1.sir hemi2.sir .. > output.rad\n", progname);
 		return(1);
 	}
 	fputs("# ", stdout);			/* copy our command */
 	printargs(argc, argv, stdout);
 						/* evaluate BSDF */
 	if (inpXML) {
-		SDclearBSDF(&myBSDF, argv[1]);
-		if (SDreportError(SDloadFile(&myBSDF, argv[1]), stderr))
+		SDclearBSDF(&myBSDF, argv[a]);
+		if (SDreportError(SDloadFile(&myBSDF, argv[a]), stderr))
 			return(1);
 		if (myBSDF.rf != NULL) front_comp |= SDsampR;
 		if (myBSDF.tf != NULL) front_comp |= SDsampT;
@@ -729,16 +742,21 @@ main(int argc, char *argv[])
 		if (myBSDF.tb != NULL) back_comp |= SDsampT;
 		if (!front_comp & !back_comp) {
 			fprintf(stderr, "%s: nothing to plot in '%s'\n",
-					progname, argv[1]);
+					progname, argv[a]);
 			return(1);
 		}
-		if (front_comp & SDsampR && myBSDF.rLambFront.cieY < overall_min*PI)
-			overall_min = myBSDF.rLambFront.cieY/PI;
-		if (back_comp & SDsampR && myBSDF.rLambBack.cieY < overall_min*PI)
-			overall_min = myBSDF.rLambBack.cieY/PI;
-		if ((front_comp|back_comp) & SDsampT &&
-				myBSDF.tLamb.cieY < overall_min*PI)
-			overall_min = myBSDF.tLamb.cieY/PI;
+		if (myLim[0] >= 0)
+			overall_min = myLim[0];
+		else {
+			if (front_comp & SDsampR && myBSDF.rLambFront.cieY < overall_min*PI)
+				overall_min = myBSDF.rLambFront.cieY/PI;
+			if (back_comp & SDsampR && myBSDF.rLambBack.cieY < overall_min*PI)
+				overall_min = myBSDF.rLambBack.cieY/PI;
+			if (front_comp & SDsampT && myBSDF.tLambFront.cieY < overall_min*PI)
+				overall_min = myBSDF.tLambFront.cieY/PI;
+			if (back_comp & SDsampT && myBSDF.tLambBack.cieY < overall_min*PI)
+				overall_min = myBSDF.tLambBack.cieY/PI;
+		}
 		set_minlog();
 		if (!build_wBSDF(&myBSDF))
 			return(1);
@@ -747,40 +765,51 @@ main(int argc, char *argv[])
 		else
 			strcpy(bsdf_name, myBSDF.name);
 		strcpy(bsdf_manuf, myBSDF.makr);
-		put_matBSDF(argv[1]);
+		put_matBSDF(argv[a]);
 	} else {
-		FILE	*fp;
-		for (n = 1; n < argc; n++) {
-			fp = fopen(argv[n], "rb");
-			if (fp == NULL) {
+		FILE	*fp[4];
+		if (argc > a+4) {
+			fprintf(stderr, "%s: more than 4 hemispheres!\n", progname);
+			return(1);
+		}
+		for (n = a; n < argc; n++) {
+			fp[n-a] = fopen(argv[n], "rb");
+			if (fp[n-a] == NULL) {
 				fprintf(stderr, "%s: cannot open BSDF interpolant '%s'\n",
 						progname, argv[n]);
 				return(1);
 			}
-			if (getheader(fp, rbf_headline, NULL) < 0) {
+			if (getheader(fp[n-a], rbf_headline, NULL) < 0) {
 				fprintf(stderr, "%s: bad BSDF interpolant '%s'\n",
 						progname, argv[n]);
 				return(1);
 			}
-			fclose(fp);
 		}
+		if (myLim[0] >= 0)
+			overall_min = myLim[0];
 		set_minlog();
-		for (n = 1; n < argc; n++) {
-			fp = fopen(argv[n], "rb");
-			if (!load_bsdf_rep(fp))
+		for (n = a; n < argc; n++) {
+			if (fseek(fp[n-a], 0L, SEEK_SET) < 0) {
+				fprintf(stderr, "%s: cannot seek on '%s'\n",
+						progname, argv[n]);
 				return(1);
-			fclose(fp);
+			}
+			if (!load_bsdf_rep(fp[n-a]))
+				return(1);
+			fclose(fp[n-a]);
 			if (!build_wRBF())
 				return(1);
 		}
 		put_matBSDF(NULL);
 	}
+	if (myLim[1] > myLim[0])	/* override maximum BSDF? */
+		overall_max = myLim[1];
 	put_source();			/* before hemispheres & labels */
 	put_hemispheres();
 	put_scale();
 	if (inpXML && myBSDF.mgf)
 		convert_mgf(myBSDF.mgf);
-	if (!put_BSDFs())
+	if (!put_BSDFs())		/* most of the output happens here */
 		return(1);
 	cleanup_tmp();
 	return(0);