--- ray/src/util/genBSDF.pl	2011/02/21 22:48:51	2.9
+++ ray/src/util/genBSDF.pl	2011/04/16 01:13:22	2.14
@@ -1,18 +1,20 @@
 #!/usr/bin/perl -w
-# RCSid $Id: genBSDF.pl,v 2.9 2011/02/21 22:48:51 greg Exp $
+# RCSid $Id: genBSDF.pl,v 2.14 2011/04/16 01:13:22 greg Exp $
 #
 # Compute BSDF based on geometry and material description
 #
 #	G. Ward
 #
 use strict;
+use File::Temp qw/ :mktemp  /;
 sub userror {
-	print STDERR "Usage: genBSDF [-n Nproc][-c Nsamp][-dim xmin xmax ymin ymax zmin zmax][{+|-}f][{+|-}b][{+|-}mgf][{+|-}geom] [input ..]\n";
+	print STDERR "Usage: genBSDF [-n Nproc][-c Nsamp][-r \"ropts\"][-dim xmin xmax ymin ymax zmin zmax][{+|-}f][{+|-}b][{+|-}mgf][{+|-}geom] [input ..]\n";
 	exit 1;
 }
-my $td = `mktemp -d /tmp/genBSDF.XXXXXX`;
+my $td = mkdtemp("/tmp/genBSDF.XXXXXX");
 chomp $td;
 my $nsamp = 1000;
+my $rtargs = "-w -ab 5 -ad 700 -lw 3e-6";
 my $mgfin = 0;
 my $geout = 1;
 my $nproc = 1;
@@ -23,6 +25,9 @@ my @dim;
 while ($#ARGV >= 0) {
 	if ("$ARGV[0]" =~ /^[-+]m/) {
 		$mgfin = ("$ARGV[0]" =~ /^\+/);
+	} elsif ("$ARGV[0]" eq "-r") {
+		$rtargs = "$rtargs $ARGV[1]";
+		shift @ARGV;
 	} elsif ("$ARGV[0]" =~ /^[-+]g/) {
 		$geout = ("$ARGV[0]" =~ /^\+/);
 	} elsif ("$ARGV[0]" =~ /^[-+]f/) {
@@ -57,7 +62,7 @@ if ( $mgfin ) {
 	die "Could not load MGF input\n" if ( $? );
 	system "mgf2rad $mgfscn > $radscn";
 } else {
-	system "cat @ARGV | xform -e > $radscn";
+	system "xform -e @ARGV > $radscn";
 	die "Could not load Radiance input\n" if ( $? );
 	system "rad2mgf $radscn > $mgfscn" if ( $geout );
 }
@@ -66,21 +71,13 @@ if ($#dim != 5) {
 }
 print STDERR "Warning: Device extends into room\n" if ($dim[5] > 1e-5);
 # Add receiver surfaces (rectangular)
-my $bmodnm="receiver_behind";
 my $fmodnm="receiver_face";
+my $bmodnm="receiver_behind";
 open(RADSCN, ">> $radscn");
-print RADSCN "void glow $fmodnm\n0\n0\n4 0 0 0 0\n\n";
-print RADSCN "$fmodnm polygon f_receiver\n0\n0\n12\n";
-print RADSCN "\t",$dim[0],"\t",$dim[2],"\t",$dim[5]+2e-5,"\n";
-print RADSCN "\t",$dim[0],"\t",$dim[3],"\t",$dim[5]+2e-5,"\n";
-print RADSCN "\t",$dim[1],"\t",$dim[3],"\t",$dim[5]+2e-5,"\n";
-print RADSCN "\t",$dim[1],"\t",$dim[2],"\t",$dim[5]+2e-5,"\n";
-print RADSCN "void glow $bmodnm\n0\n0\n4 0 0 0 0\n\n";
-print RADSCN "$bmodnm polygon b_receiver\n0\n0\n12\n";
-print RADSCN "\t",$dim[1],"\t",$dim[2],"\t",$dim[4]-2e-5,"\n";
-print RADSCN "\t",$dim[1],"\t",$dim[3],"\t",$dim[4]-2e-5,"\n";
-print RADSCN "\t",$dim[0],"\t",$dim[3],"\t",$dim[4]-2e-5,"\n";
-print RADSCN "\t",$dim[0],"\t",$dim[2],"\t",$dim[4]-2e-5,"\n";
+print RADSCN "void glow $fmodnm\n0\n0\n4 1 1 1 0\n\n";
+print RADSCN "$fmodnm source f_receiver\n0\n0\n4 0 0 1 180\n";
+print RADSCN "void glow $bmodnm\n0\n0\n4 1 1 1 0\n\n";
+print RADSCN "$bmodnm source b_receiver\n0\n0\n4 0 0 -1 180\n";
 close RADSCN;
 # Generate octree
 system "oconv -w $radscn > $octree";
@@ -110,6 +107,7 @@ KprojOmega = PI * if(Kbin-.5,
 # Compute Klems bin from exiting ray direction (forward or backward)
 my $kcal = '
 DEGREE : PI/180;
+abs(x) : if(x, x, -x);
 Acos(x) : 1/DEGREE * if(x-1, 0, if(-1-x, 0, acos(x)));
 posangle(a) : if(-a, a + 2*PI, a);
 Atan2(y,x) : 1/DEGREE * posangle(atan2(y,x));
@@ -130,7 +128,7 @@ kbin2(pol,azi) = select(kfindrow(1, pol),
 		kaccum(7) + kazn(azi,360/knaz(8)),
 		kaccum(8) + kazn(azi,360/knaz(9))
 	);
-kbin = if(Dz, kbin2(Acos(Dz),Atan2(Dy,Dx)), kbin2(Acos(-Dz),Atan2(-Dy,-Dx)));
+kbin = kbin2(Acos(abs(Dz)),Atan2(Dy,Dx));
 ';
 my $ndiv = 145;
 my $nx = int(sqrt($nsamp*($dim[1]-$dim[0])/($dim[3]-$dim[2])) + .5);
@@ -142,7 +140,6 @@ my @rfarr;
 my @tbarr;
 my @rbarr;
 my $cmd;
-my $rtargs = "-w -ab 5 -ad 700 -lw 3e-6";
 my $rtcmd = "rtcontrib -h -ff -fo -n $nproc -c $nsamp " .
 	"-e '$kcal' -b kbin -bn $ndiv " .
 	"-o '$td/%s.flt' -m $fmodnm -m $bmodnm $rtargs $octree";
@@ -151,11 +148,11 @@ my $rccmd = "rcalc -e '$tcal' " .
 	q{-if3 -e '$1=(0.265*$1+0.670*$2+0.065*$3)/KprojOmega'};
 if ( $doforw ) {
 $cmd = "cnt $ndiv $ny $nx | rcalc -of -e '$tcal' " .
-	"-e 'xp=(\$3+rand(.35*recno-15))*(($dim[1]-$dim[0])/$nx)+$dim[0]' " .
-	"-e 'yp=(\$2+rand(.86*recno+11))*(($dim[3]-$dim[2])/$ny)+$dim[2]' " .
-	"-e 'zp:$dim[4]-1e-5' " .
-	q{-e 'Kbin=$1;x1=rand(1.21*recno+2.75);x2=rand(-3.55*recno-7.57)' } .
-	q{-e '$1=xp;$2=yp;$3=zp;$4=Dx;$5=Dy;$6=Dz' } .
+	"-e 'xp=(\$3+rand(.12*recno+288))*(($dim[1]-$dim[0])/$nx)+$dim[0]' " .
+	"-e 'yp=(\$2+rand(.37*recno-44))*(($dim[3]-$dim[2])/$ny)+$dim[2]' " .
+	"-e 'zp:$dim[4]' " .
+	q{-e 'Kbin=$1;x1=rand(2.75*recno+3.1);x2=rand(-2.01*recno-3.37)' } .
+	q{-e '$1=xp-Dx;$2=yp-Dy;$3=zp-Dz;$4=Dx;$5=Dy;$6=Dz' } .
 	"| $rtcmd";
 system "$cmd" || die "Failure running: $cmd\n";
 @tfarr = `$rccmd $td/$fmodnm.flt`;
@@ -167,9 +164,9 @@ if ( $doback ) {
 $cmd = "cnt $ndiv $ny $nx | rcalc -of -e '$tcal' " .
 	"-e 'xp=(\$3+rand(.35*recno-15))*(($dim[1]-$dim[0])/$nx)+$dim[0]' " .
 	"-e 'yp=(\$2+rand(.86*recno+11))*(($dim[3]-$dim[2])/$ny)+$dim[2]' " .
-	"-e 'zp:$dim[5]+1e-5' " .
+	"-e 'zp:$dim[5]' " .
 	q{-e 'Kbin=$1;x1=rand(1.21*recno+2.75);x2=rand(-3.55*recno-7.57)' } .
-	q{-e '$1=xp;$2=yp;$3=zp;$4=-Dx;$5=-Dy;$6=-Dz' } .
+	q{-e '$1=xp-Dx;$2=yp-Dy;$3=zp+Dz;$4=Dx;$5=Dy;$6=-Dz' } .
 	"| $rtcmd";
 system "$cmd" || die "Failure running: $cmd\n";
 @tbarr = `$rccmd $td/$bmodnm.flt`;
@@ -196,7 +193,9 @@ print "\t\t\t<DeviceType>Integral</DeviceType>\n";
 if ( $geout ) {
 	print "\t\t\t<Geometry format=\"MGF\" unit=\"Meter\">\n";
 	printf "xf -t %.6f %.6f 0\n", -($dim[0]+$dim[1])/2, -($dim[2]+$dim[3])/2;
-	system "cat $mgfscn";
+	open(MGFSCN, "< $mgfscn");
+	while (<MGFSCN>) { print $_; }
+	close MGFSCN;
 	print "xf\n";
 	print "\t\t\t</Geometry>\n";
 }