| 1 | greg | 2.1 | #!/usr/bin/perl -w | 
| 2 | greg | 2.89 | # RCSid $Id: genBSDF.pl,v 2.88 2022/01/24 17:49:32 greg Exp $ | 
| 3 | greg | 2.1 | # | 
| 4 |  |  | # Compute BSDF based on geometry and material description | 
| 5 |  |  | # | 
| 6 |  |  | #       G. Ward | 
| 7 |  |  | # | 
| 8 |  |  | use strict; | 
| 9 | greg | 2.47 | my $windoz = ($^O eq "MSWin32" or $^O eq "MSWin64"); | 
| 10 | greg | 2.13 | use File::Temp qw/ :mktemp  /; | 
| 11 | greg | 2.1 | sub userror { | 
| 12 | greg | 2.76 | print STDERR "Usage: genBSDF [-n Nproc][-c Nsamp][-W][-t{3|4} Nlog2][-r \"ropts\"][-s \"x=string;y=string\"][-dim xmin xmax ymin ymax zmin zmax][{+|-}C][{+|-}a][{+|-}f][{+|-}b][{+|-}mgf][{+|-}geom units] [input ..]\n"; | 
| 13 | greg | 2.1 | exit 1; | 
| 14 |  |  | } | 
| 15 | greg | 2.49 | my ($td,$radscn,$mgfscn,$octree,$fsender,$bsender,$receivers,$facedat,$behinddat,$rmtmp); | 
| 16 | greg | 2.68 | my ($tf,$rf,$tb,$rb,$tfx,$rfx,$tbx,$rbx,$tfz,$rfz,$tbz,$rbz,$cph); | 
| 17 | greg | 2.66 | my ($curphase, $recovery); | 
| 18 |  |  | if ($#ARGV == 1 && "$ARGV[0]" =~ /^-rec/) { | 
| 19 |  |  | $td = $ARGV[1]; | 
| 20 |  |  | open(MYAVH, "< $td/savedARGV.txt") or die "$td: invalid path\n"; | 
| 21 |  |  | @ARGV = <MYAVH>; | 
| 22 |  |  | close MYAVH; | 
| 23 |  |  | chomp @ARGV; | 
| 24 | greg | 2.69 | $recovery = 0; | 
| 25 | greg | 2.66 | if (open(MYPH, "< $td/phase.txt")) { | 
| 26 |  |  | while (<MYPH>) { | 
| 27 |  |  | chomp($recovery = $_); | 
| 28 |  |  | } | 
| 29 |  |  | close MYPH; | 
| 30 |  |  | } | 
| 31 |  |  | } elsif ($windoz) { | 
| 32 | greg | 2.47 | my $tmploc = `echo \%TMP\%`; | 
| 33 | greg | 2.66 | chomp $tmploc; | 
| 34 | greg | 2.47 | $td = mkdtemp("$tmploc\\genBSDF.XXXXXX"); | 
| 35 | greg | 2.66 | } else { | 
| 36 |  |  | $td = mkdtemp("/tmp/genBSDF.XXXXXX"); | 
| 37 |  |  | chomp $td; | 
| 38 |  |  | } | 
| 39 |  |  | if ($windoz) { | 
| 40 | greg | 2.47 | $radscn = "$td\\device.rad"; | 
| 41 |  |  | $mgfscn = "$td\\device.mgf"; | 
| 42 |  |  | $octree = "$td\\device.oct"; | 
| 43 | greg | 2.49 | $fsender = "$td\\fsender.rad"; | 
| 44 |  |  | $bsender = "$td\\bsender.rad"; | 
| 45 |  |  | $receivers = "$td\\receivers.rad"; | 
| 46 |  |  | $facedat = "$td\\face.dat"; | 
| 47 |  |  | $behinddat = "$td\\behind.dat"; | 
| 48 |  |  | $tf = "$td\\tf.dat"; | 
| 49 |  |  | $rf = "$td\\rf.dat"; | 
| 50 |  |  | $tb = "$td\\tb.dat"; | 
| 51 |  |  | $rb = "$td\\rb.dat"; | 
| 52 |  |  | $tfx = "$td\\tfx.dat"; | 
| 53 |  |  | $rfx = "$td\\rfx.dat"; | 
| 54 |  |  | $tbx = "$td\\tbx.dat"; | 
| 55 |  |  | $rbx = "$td\\rbx.dat"; | 
| 56 |  |  | $tfz = "$td\\tfz.dat"; | 
| 57 |  |  | $rfz = "$td\\rfz.dat"; | 
| 58 |  |  | $tbz = "$td\\tbz.dat"; | 
| 59 |  |  | $rbz = "$td\\rbz.dat"; | 
| 60 | greg | 2.66 | $cph = "$td\\phase.txt"; | 
| 61 | greg | 2.63 | $rmtmp = "rd /S /Q $td"; | 
| 62 | greg | 2.49 | } else { | 
| 63 | greg | 2.47 | $radscn = "$td/device.rad"; | 
| 64 |  |  | $mgfscn = "$td/device.mgf"; | 
| 65 |  |  | $octree = "$td/device.oct"; | 
| 66 | greg | 2.49 | $fsender = "$td/fsender.rad"; | 
| 67 |  |  | $bsender = "$td/bsender.rad"; | 
| 68 |  |  | $receivers = "$td/receivers.rad"; | 
| 69 |  |  | $facedat = "$td/face.dat"; | 
| 70 |  |  | $behinddat = "$td/behind.dat"; | 
| 71 |  |  | $tf = "$td/tf.dat"; | 
| 72 |  |  | $rf = "$td/rf.dat"; | 
| 73 |  |  | $tb = "$td/tb.dat"; | 
| 74 |  |  | $rb = "$td/rb.dat"; | 
| 75 |  |  | $tfx = "$td/tfx.dat"; | 
| 76 |  |  | $rfx = "$td/rfx.dat"; | 
| 77 |  |  | $tbx = "$td/tbx.dat"; | 
| 78 |  |  | $rbx = "$td/rbx.dat"; | 
| 79 |  |  | $tfz = "$td/tfz.dat"; | 
| 80 |  |  | $rfz = "$td/rfz.dat"; | 
| 81 |  |  | $tbz = "$td/tbz.dat"; | 
| 82 |  |  | $rbz = "$td/rbz.dat"; | 
| 83 | greg | 2.66 | $cph = "$td/phase.txt"; | 
| 84 | greg | 2.47 | $rmtmp = "rm -rf $td"; | 
| 85 |  |  | } | 
| 86 | greg | 2.20 | my @savedARGV = @ARGV; | 
| 87 | greg | 2.74 | my $rfluxmtx = "rfluxmtx -ab 5 -ad 700 -lw 3e-6 -w-"; | 
| 88 | greg | 2.49 | my $wrapper = "wrapBSDF"; | 
| 89 | greg | 2.15 | my $tensortree = 0; | 
| 90 |  |  | my $ttlog2 = 4; | 
| 91 | greg | 2.76 | my $dorecip = 1; | 
| 92 | greg | 2.25 | my $nsamp = 2000; | 
| 93 | greg | 2.1 | my $mgfin = 0; | 
| 94 |  |  | my $geout = 1; | 
| 95 |  |  | my $nproc = 1; | 
| 96 | greg | 2.49 | my $docolor = 0; | 
| 97 | greg | 2.9 | my $doforw = 0; | 
| 98 |  |  | my $doback = 1; | 
| 99 | greg | 2.33 | my $pctcull = 90; | 
| 100 | greg | 2.41 | my $gunit = "meter"; | 
| 101 | greg | 2.49 | my $curspec = "Visible"; | 
| 102 | greg | 2.1 | my @dim; | 
| 103 |  |  | # Get options | 
| 104 |  |  | while ($#ARGV >= 0) { | 
| 105 |  |  | if ("$ARGV[0]" =~ /^[-+]m/) { | 
| 106 |  |  | $mgfin = ("$ARGV[0]" =~ /^\+/); | 
| 107 | greg | 2.10 | } elsif ("$ARGV[0]" eq "-r") { | 
| 108 | greg | 2.49 | $rfluxmtx .= " $ARGV[1]"; | 
| 109 | greg | 2.10 | shift @ARGV; | 
| 110 | greg | 2.1 | } elsif ("$ARGV[0]" =~ /^[-+]g/) { | 
| 111 |  |  | $geout = ("$ARGV[0]" =~ /^\+/); | 
| 112 | greg | 2.22 | $gunit = $ARGV[1]; | 
| 113 |  |  | if ($gunit !~ /^(?i)(meter|foot|inch|centimeter|millimeter)$/) { | 
| 114 |  |  | die "Illegal geometry unit '$gunit': must be meter, foot, inch, centimeter, or millimeter\n"; | 
| 115 |  |  | } | 
| 116 |  |  | shift @ARGV; | 
| 117 | greg | 2.49 | } elsif ("$ARGV[0]" =~ /^[-+]C/) { | 
| 118 |  |  | $docolor = ("$ARGV[0]" =~ /^\+/); | 
| 119 | greg | 2.78 | } elsif ("$ARGV[0]" =~ /^[-+]a/) { | 
| 120 | greg | 2.76 | $dorecip = ("$ARGV[0]" =~ /^\+/); | 
| 121 | greg | 2.9 | } elsif ("$ARGV[0]" =~ /^[-+]f/) { | 
| 122 |  |  | $doforw = ("$ARGV[0]" =~ /^\+/); | 
| 123 |  |  | } elsif ("$ARGV[0]" =~ /^[-+]b/) { | 
| 124 |  |  | $doback = ("$ARGV[0]" =~ /^\+/); | 
| 125 | greg | 2.27 | } elsif ("$ARGV[0]" eq "-t") { | 
| 126 | greg | 2.33 | # Use value < 0 for rttree_reduce bypass | 
| 127 | greg | 2.27 | $pctcull = $ARGV[1]; | 
| 128 | greg | 2.55 | if ($pctcull >= 100) { | 
| 129 |  |  | die "Illegal -t culling percentage, must be < 100\n"; | 
| 130 |  |  | } | 
| 131 | greg | 2.27 | shift @ARGV; | 
| 132 | greg | 2.15 | } elsif ("$ARGV[0]" =~ /^-t[34]$/) { | 
| 133 |  |  | $tensortree = substr($ARGV[0], 2, 1); | 
| 134 |  |  | $ttlog2 = $ARGV[1]; | 
| 135 |  |  | shift @ARGV; | 
| 136 | greg | 2.55 | } elsif ("$ARGV[0]" eq "-s") { | 
| 137 | greg | 2.49 | $wrapper .= " -f \"$ARGV[1]\""; | 
| 138 |  |  | shift @ARGV; | 
| 139 |  |  | } elsif ("$ARGV[0]" eq "-W") { | 
| 140 |  |  | $wrapper .= " -W"; | 
| 141 | greg | 2.1 | } elsif ("$ARGV[0]" eq "-c") { | 
| 142 |  |  | $nsamp = $ARGV[1]; | 
| 143 |  |  | shift @ARGV; | 
| 144 |  |  | } elsif ("$ARGV[0]" eq "-n") { | 
| 145 |  |  | $nproc = $ARGV[1]; | 
| 146 |  |  | shift @ARGV; | 
| 147 |  |  | } elsif ("$ARGV[0]" =~ /^-d/) { | 
| 148 |  |  | userror() if ($#ARGV < 6); | 
| 149 | greg | 2.8 | @dim = @ARGV[1..6]; | 
| 150 | greg | 2.1 | shift @ARGV for (1..6); | 
| 151 |  |  | } elsif ("$ARGV[0]" =~ /^[-+]./) { | 
| 152 |  |  | userror(); | 
| 153 |  |  | } else { | 
| 154 |  |  | last; | 
| 155 |  |  | } | 
| 156 |  |  | shift @ARGV; | 
| 157 |  |  | } | 
| 158 | greg | 2.9 | # Check that we're actually being asked to do something | 
| 159 | greg | 2.15 | die "Must have at least one of +forward or +backward\n" if (!$doforw && !$doback); | 
| 160 | greg | 2.49 | $wrapper .= $tensortree ? " -a t$tensortree" : " -a kf -c"; | 
| 161 |  |  | $wrapper .= " -u $gunit"; | 
| 162 | greg | 2.67 | if ( !defined $recovery ) { | 
| 163 | greg | 2.66 | # Issue warning for unhandled reciprocity case | 
| 164 |  |  | print STDERR "Warning: recommend both +forward and +backward with -t3\n" if | 
| 165 |  |  | ($tensortree==3 && !($doforw && $doback)); | 
| 166 | greg | 2.67 | # Get scene description | 
| 167 | greg | 2.66 | if ( $mgfin ) { | 
| 168 |  |  | system "mgf2rad @ARGV > $radscn"; | 
| 169 |  |  | die "Could not load MGF input\n" if ( $? ); | 
| 170 |  |  | } else { | 
| 171 |  |  | system "xform -e @ARGV > $radscn"; | 
| 172 |  |  | die "Could not load Radiance input\n" if ( $? ); | 
| 173 |  |  | } | 
| 174 | greg | 2.1 | } | 
| 175 | greg | 2.67 | if ( $#dim != 5 ) { | 
| 176 | greg | 2.74 | @dim = split ' ', `getbbox -h -w $radscn`; | 
| 177 | greg | 2.1 | } | 
| 178 | greg | 2.67 | die "Device entirely inside room!\n" if ( $dim[4] >= 0 ); | 
| 179 |  |  | if ( $dim[5] > 1e-5 ) { | 
| 180 | greg | 2.64 | print STDERR "Warning: Device extends into room\n"; | 
| 181 | greg | 2.67 | } elsif ( $dim[5]*$dim[5] > .01*($dim[1]-$dim[0])*($dim[3]-$dim[2]) ) { | 
| 182 | greg | 2.64 | print STDERR "Warning: Device far behind Z==0 plane\n"; | 
| 183 |  |  | } | 
| 184 |  |  | # Assume Zmax==0 to derive thickness so pkgBSDF will work | 
| 185 |  |  | $wrapper .= ' -f "t=' . (-$dim[4]) . ';w=' . ($dim[1] - $dim[0]) . | 
| 186 | greg | 2.49 | ';h=' . ($dim[3] - $dim[2]) . '"'; | 
| 187 | greg | 2.67 | $wrapper .= " -g $mgfscn" if ( $geout ); | 
| 188 | greg | 2.66 | # Calculate CIE (u',v') from Radiance RGB: | 
| 189 |  |  | my $CIEuv =     'Xi=.5141*Ri+.3239*Gi+.1620*Bi;' . | 
| 190 |  |  | 'Yi=.2651*Ri+.6701*Gi+.0648*Bi;' . | 
| 191 |  |  | 'Zi=.0241*Ri+.1229*Gi+.8530*Bi;' . | 
| 192 | greg | 2.80 | 'den=Xi+15*Yi+3*Zi;' . | 
| 193 |  |  | 'uprime=if(Yi,4*Xi/den,4/19);' . | 
| 194 |  |  | 'vprime=if(Yi,9*Yi/den,9/19);' ; | 
| 195 | greg | 2.56 | my $FEPS = 1e-5; | 
| 196 | greg | 2.66 | my $ns = 2**$ttlog2; | 
| 197 | greg | 2.46 | my $nx = int(sqrt($nsamp*($dim[1]-$dim[0])/($dim[3]-$dim[2])) + 1); | 
| 198 |  |  | my $ny = int($nsamp/$nx + 1); | 
| 199 | greg | 2.15 | $nsamp = $nx * $ny; | 
| 200 | greg | 2.49 | $rfluxmtx .= " -n $nproc -c $nsamp"; | 
| 201 | greg | 2.67 | if ( !defined $recovery ) { | 
| 202 | greg | 2.66 | open(MYAVH, "> $td/savedARGV.txt"); | 
| 203 |  |  | foreach (@savedARGV) { | 
| 204 |  |  | print MYAVH "$_\n"; | 
| 205 |  |  | } | 
| 206 |  |  | close MYAVH; | 
| 207 |  |  | # Generate octree | 
| 208 | greg | 2.81 | system "oconv -w $radscn > $octree"; | 
| 209 | greg | 2.66 | die "Could not compile scene\n" if ( $? ); | 
| 210 |  |  | # Add MGF description if requested | 
| 211 |  |  | if ( $geout ) { | 
| 212 |  |  | open(MGFSCN, "> $mgfscn"); | 
| 213 |  |  | printf MGFSCN "xf -t %.6f %.6f 0\n", -($dim[0]+$dim[1])/2, -($dim[2]+$dim[3])/2; | 
| 214 |  |  | close MGFSCN; | 
| 215 |  |  | if ( $mgfin ) { | 
| 216 | greg | 2.89 | die "+mgf requires input file with +geom\n" if ($#ARGV < 0); | 
| 217 | greg | 2.66 | system qq{mgfilt "#,o,xf,c,cxy,cspec,cmix,m,sides,rd,td,rs,ts,ir,v,p,n,f,fh,sph,cyl,cone,prism,ring,torus" @ARGV >> $mgfscn}; | 
| 218 |  |  | } else { | 
| 219 |  |  | system "rad2mgf $radscn >> $mgfscn"; | 
| 220 |  |  | } | 
| 221 |  |  | open(MGFSCN, ">> $mgfscn"); | 
| 222 |  |  | print MGFSCN "xf\n"; | 
| 223 |  |  | close MGFSCN; | 
| 224 |  |  | } | 
| 225 |  |  | # Create receiver & sender surfaces (rectangular) | 
| 226 |  |  | open(RADSCN, "> $receivers"); | 
| 227 |  |  | print RADSCN '#@rfluxmtx ' . ($tensortree ? "h=-sc$ns\n" : "h=-kf\n"); | 
| 228 |  |  | print RADSCN '#@rfluxmtx ' . "u=-Y o=$facedat\n\n"; | 
| 229 |  |  | print RADSCN "void glow receiver_face\n0\n0\n4 1 1 1 0\n\n"; | 
| 230 |  |  | print RADSCN "receiver_face source f_receiver\n0\n0\n4 0 0 1 180\n\n"; | 
| 231 |  |  | print RADSCN '#@rfluxmtx ' . ($tensortree ? "h=+sc$ns\n" : "h=+kf\n"); | 
| 232 |  |  | print RADSCN '#@rfluxmtx ' . "u=-Y o=$behinddat\n\n"; | 
| 233 |  |  | print RADSCN "void glow receiver_behind\n0\n0\n4 1 1 1 0\n\n"; | 
| 234 |  |  | print RADSCN "receiver_behind source b_receiver\n0\n0\n4 0 0 -1 180\n"; | 
| 235 | greg | 2.49 | close RADSCN; | 
| 236 | greg | 2.66 | # Prepare sender surfaces | 
| 237 |  |  | if ( $tensortree != 3 ) {       # Isotropic tensor tree is exception | 
| 238 |  |  | open (RADSCN, "> $fsender"); | 
| 239 |  |  | print RADSCN '#@rfluxmtx u=-Y ' . ($tensortree ? "h=-sc$ns\n\n" : "h=-kf\n\n"); | 
| 240 |  |  | print RADSCN "void polygon fwd_sender\n0\n0\n12\n"; | 
| 241 |  |  | printf RADSCN "\t%e\t%e\t%e\n", $dim[0], $dim[2], $dim[4]-$FEPS; | 
| 242 |  |  | printf RADSCN "\t%e\t%e\t%e\n", $dim[0], $dim[3], $dim[4]-$FEPS; | 
| 243 |  |  | printf RADSCN "\t%e\t%e\t%e\n", $dim[1], $dim[3], $dim[4]-$FEPS; | 
| 244 |  |  | printf RADSCN "\t%e\t%e\t%e\n", $dim[1], $dim[2], $dim[4]-$FEPS; | 
| 245 |  |  | close RADSCN; | 
| 246 |  |  | open (RADSCN, "> $bsender"); | 
| 247 |  |  | print RADSCN '#@rfluxmtx u=-Y ' . ($tensortree ? "h=+sc$ns\n\n" : "h=+kf\n\n"); | 
| 248 |  |  | print RADSCN "void polygon bwd_sender\n0\n0\n12\n"; | 
| 249 |  |  | printf RADSCN "\t%e\t%e\t%e\n", $dim[0], $dim[2], $dim[5]+$FEPS; | 
| 250 |  |  | printf RADSCN "\t%e\t%e\t%e\n", $dim[1], $dim[2], $dim[5]+$FEPS; | 
| 251 |  |  | printf RADSCN "\t%e\t%e\t%e\n", $dim[1], $dim[3], $dim[5]+$FEPS; | 
| 252 |  |  | printf RADSCN "\t%e\t%e\t%e\n", $dim[0], $dim[3], $dim[5]+$FEPS; | 
| 253 |  |  | close RADSCN; | 
| 254 |  |  | } | 
| 255 |  |  | print STDERR "Recover using: $0 -recover $td\n"; | 
| 256 |  |  | } | 
| 257 |  |  | $curphase = 0; | 
| 258 | greg | 2.15 | # Create data segments (all the work happens here) | 
| 259 |  |  | if ( $tensortree ) { | 
| 260 |  |  | do_tree_bsdf(); | 
| 261 |  |  | } else { | 
| 262 |  |  | do_matrix_bsdf(); | 
| 263 |  |  | } | 
| 264 | greg | 2.49 | # Output XML | 
| 265 | greg | 2.67 | # print STDERR "Running: $wrapper\n"; | 
| 266 | greg | 2.62 | system "$wrapper -C \"Created by: genBSDF @savedARGV\""; | 
| 267 | greg | 2.49 | die "Could not wrap BSDF data\n" if ( $? ); | 
| 268 | greg | 2.15 | # Clean up temporary files and exit | 
| 269 | greg | 2.53 | exec $rmtmp; | 
| 270 | greg | 2.16 |  | 
| 271 | greg | 2.67 | #============== End of main program segment ==============# | 
| 272 |  |  |  | 
| 273 |  |  | # Function to determine if next phase should be skipped or recovered | 
| 274 |  |  | sub do_phase { | 
| 275 |  |  | $curphase++; | 
| 276 | greg | 2.70 | if ( defined $recovery ) { | 
| 277 |  |  | if ( $recovery > $curphase ) { return 0; } | 
| 278 |  |  | if ( $recovery == $curphase ) { return -1; } | 
| 279 | greg | 2.67 | } | 
| 280 | greg | 2.73 | open(MYPH, ">> $td/phase.txt"); | 
| 281 | greg | 2.67 | print MYPH "$curphase\n"; | 
| 282 | greg | 2.73 | close MYPH; | 
| 283 | greg | 2.67 | return 1; | 
| 284 |  |  | } | 
| 285 | greg | 2.15 |  | 
| 286 | greg | 2.70 | # Check if we are in active phase (not skipping parts) | 
| 287 |  |  | sub active_phase { | 
| 288 |  |  | if ( defined $recovery ) { | 
| 289 |  |  | if ( $recovery > $curphase ) { return 0; } | 
| 290 |  |  | if ( $recovery == $curphase ) { return -1; } | 
| 291 |  |  | } | 
| 292 |  |  | return 1; | 
| 293 |  |  | } | 
| 294 |  |  |  | 
| 295 |  |  | # Function to run program and check output if in active phase | 
| 296 |  |  | sub run_check { | 
| 297 |  |  | if ( !active_phase() ) { return; } | 
| 298 |  |  | my $cmd = shift; | 
| 299 |  |  | # print STDERR "Running: $cmd\n"; | 
| 300 |  |  | system $cmd; | 
| 301 |  |  | die "Failure running: $cmd\n" if ( $? ); | 
| 302 |  |  | } | 
| 303 |  |  |  | 
| 304 | greg | 2.15 | #++++++++++++++ Tensor tree BSDF generation ++++++++++++++# | 
| 305 |  |  | sub do_tree_bsdf { | 
| 306 | greg | 2.49 |  | 
| 307 |  |  | # Run rfluxmtx processes to compute each side | 
| 308 |  |  | do_ttree_dir(0) if ( $doback ); | 
| 309 |  |  | do_ttree_dir(1) if ( $doforw ); | 
| 310 | greg | 2.36 |  | 
| 311 | greg | 2.15 | }       # end of sub do_tree_bsdf() | 
| 312 |  |  |  | 
| 313 | greg | 2.49 | # Call rfluxmtx and process tensor tree BSDF for the given direction | 
| 314 |  |  | sub do_ttree_dir { | 
| 315 | greg | 2.15 | my $forw = shift; | 
| 316 | greg | 2.70 | my $dop = do_phase(); | 
| 317 |  |  | my $r = ($dop < 0) ? " -r" : ""; | 
| 318 | greg | 2.47 | my $cmd; | 
| 319 | greg | 2.15 | if ( $tensortree == 3 ) { | 
| 320 |  |  | # Isotropic BSDF | 
| 321 | greg | 2.36 | my $ns2 = $ns / 2; | 
| 322 | greg | 2.47 | if ($windoz) { | 
| 323 |  |  | $cmd = "cnt $ns2 $ny $nx " . | 
| 324 |  |  | qq{| rcalc -e "r1=rand(.8681*recno-.673892)" } . | 
| 325 |  |  | qq{-e "r2=rand(-5.37138*recno+67.1737811)" } . | 
| 326 |  |  | qq{-e "r3=rand(+3.17603772*recno+83.766771)" } . | 
| 327 | greg | 2.82 | qq{-e "r4=rand(-1.5839226*recno-59.82712)" } . | 
| 328 | greg | 2.84 | qq{-e "odds(n):if(.5*n-floor(.5*n)-.25,-1,1)"} . | 
| 329 |  |  | qq{-e "Dx=1-(\$1+r1)/$ns2" } . | 
| 330 | greg | 2.86 | qq{-e "Dy=min(1/$ns,sqrt(1-Dx*Dx))*odds(\$1)*r2" } . | 
| 331 | greg | 2.82 | qq{-e "Dz=sqrt(1-Dx*Dx-Dy*Dy)" } . | 
| 332 | greg | 2.47 | qq{-e "xp=(\$3+r2)*(($dim[1]-$dim[0])/$nx)+$dim[0]" } . | 
| 333 |  |  | qq{-e "yp=(\$2+r3)*(($dim[3]-$dim[2])/$ny)+$dim[2]" } . | 
| 334 |  |  | qq{-e "zp=$dim[5-$forw]" -e "myDz=Dz*($forw*2-1)" } . | 
| 335 |  |  | qq{-e "\$1=xp-Dx;\$2=yp-Dy;\$3=zp-myDz" } . | 
| 336 |  |  | qq{-e "\$4=Dx;\$5=Dy;\$6=myDz" } . | 
| 337 | greg | 2.66 | "| $rfluxmtx$r -fa -y $ns2 - $receivers -i $octree"; | 
| 338 | greg | 2.47 | } else { | 
| 339 |  |  | $cmd = "cnt $ns2 $ny $nx " . | 
| 340 | greg | 2.75 | qq{| rcalc -e 'r1=rand(.8681*recno-.673892)' } . | 
| 341 |  |  | qq{-e 'r2=rand(-5.37138*recno+67.1737811)' } . | 
| 342 |  |  | qq{-e 'r3=rand(+3.17603772*recno+83.766771)' } . | 
| 343 | greg | 2.82 | qq{-e 'r4=rand(-1.5839226*recno-59.82712)' } . | 
| 344 | greg | 2.84 | qq{-e 'odds(n):if(.5*n-floor(.5*n)-.25,-1,1)' } . | 
| 345 |  |  | qq{-e 'Dx=1-(\$1+r1)/$ns2' } . | 
| 346 | greg | 2.86 | qq{-e 'Dy=min(1/$ns,sqrt(1-Dx*Dx))*odds(\$1)*r2' } . | 
| 347 | greg | 2.82 | qq{-e 'Dz=sqrt(1-Dx*Dx-Dy*Dy)' } . | 
| 348 |  |  | qq{-e 'xp=(\$3+r3)*(($dim[1]-$dim[0])/$nx)+$dim[0]' } . | 
| 349 |  |  | qq{-e 'yp=(\$2+r4)*(($dim[3]-$dim[2])/$ny)+$dim[2]' } . | 
| 350 | greg | 2.75 | qq{-e 'zp=$dim[5-$forw]' -e 'myDz=Dz*($forw*2-1)' } . | 
| 351 | greg | 2.47 | qq{-e '\$1=xp-Dx;\$2=yp-Dy;\$3=zp-myDz' } . | 
| 352 |  |  | qq{-e '\$4=Dx;\$5=Dy;\$6=myDz' -of } . | 
| 353 | greg | 2.66 | "| $rfluxmtx$r -h -ff -y $ns2 - $receivers -i $octree"; | 
| 354 | greg | 2.47 | } | 
| 355 | greg | 2.15 | } else { | 
| 356 |  |  | # Anisotropic BSDF | 
| 357 | greg | 2.49 | my $sender = ($bsender,$fsender)[$forw]; | 
| 358 | greg | 2.47 | if ($windoz) { | 
| 359 | greg | 2.66 | $cmd = "$rfluxmtx$r -fa $sender $receivers -i $octree"; | 
| 360 | greg | 2.47 | } else { | 
| 361 | greg | 2.66 | $cmd = "$rfluxmtx$r -h -ff $sender $receivers -i $octree"; | 
| 362 | greg | 2.47 | } | 
| 363 | greg | 2.15 | } | 
| 364 | greg | 2.71 | if ( $dop ) { | 
| 365 |  |  | # print STDERR "Running: $cmd\n"; | 
| 366 |  |  | system $cmd; | 
| 367 |  |  | die "Failure running rfluxmtx" if ( $? ); | 
| 368 |  |  | } | 
| 369 | greg | 2.36 | ttree_out($forw); | 
| 370 | greg | 2.49 | }       # end of do_ttree_dir() | 
| 371 | greg | 2.15 |  | 
| 372 | greg | 2.49 | # Simplify and store tensor tree results | 
| 373 | greg | 2.15 | sub ttree_out { | 
| 374 |  |  | my $forw = shift; | 
| 375 | greg | 2.49 | my ($refldat,$transdat); | 
| 376 |  |  | if ( $forw ) { | 
| 377 |  |  | $transdat = $facedat; | 
| 378 |  |  | $refldat = $behinddat; | 
| 379 | greg | 2.47 | } else { | 
| 380 | greg | 2.49 | $transdat = $behinddat; | 
| 381 |  |  | $refldat = $facedat; | 
| 382 | greg | 2.47 | } | 
| 383 | greg | 2.49 | # Only output one transmitted anisotropic distribution, preferring backwards | 
| 384 |  |  | if ( !$forw || !$doback || $tensortree==3 ) { | 
| 385 |  |  | my $ttyp = ("tb","tf")[$forw]; | 
| 386 |  |  | ttree_comp($ttyp, "Visible", $transdat, ($tb,$tf)[$forw]); | 
| 387 |  |  | if ( $docolor ) { | 
| 388 | greg | 2.53 | ttree_comp($ttyp, "CIE-u", $transdat, ($tbx,$tfx)[$forw]); | 
| 389 |  |  | ttree_comp($ttyp, "CIE-v", $transdat, ($tbz,$tfz)[$forw]); | 
| 390 | greg | 2.49 | } | 
| 391 | greg | 2.47 | } | 
| 392 | greg | 2.49 | # Output reflection | 
| 393 |  |  | my $rtyp = ("rb","rf")[$forw]; | 
| 394 |  |  | ttree_comp($rtyp, "Visible", $refldat, ($rb,$rf)[$forw]); | 
| 395 |  |  | if ( $docolor ) { | 
| 396 | greg | 2.53 | ttree_comp($rtyp, "CIE-u", $refldat, ($rbx,$rfx)[$forw]); | 
| 397 |  |  | ttree_comp($rtyp, "CIE-v", $refldat, ($rbz,$rfz)[$forw]); | 
| 398 | greg | 2.34 | } | 
| 399 | greg | 2.49 | }       # end of ttree_out() | 
| 400 |  |  |  | 
| 401 |  |  | # Call rttree_reduce on the given component | 
| 402 |  |  | sub ttree_comp { | 
| 403 |  |  | my $typ = shift; | 
| 404 |  |  | my $spec = shift; | 
| 405 |  |  | my $src = shift; | 
| 406 |  |  | my $dest = shift; | 
| 407 |  |  | my $cmd; | 
| 408 | greg | 2.47 | if ($windoz) { | 
| 409 | greg | 2.49 | if ("$spec" eq "Visible") { | 
| 410 |  |  | $cmd = qq{rcalc -e "Omega:PI/($ns*$ns)" } . | 
| 411 | greg | 2.53 | q{-e "Ri=$1;Gi=$2;Bi=$3" } . | 
| 412 |  |  | qq{-e "$CIEuv" } . | 
| 413 |  |  | q{-e "$1=Yi/Omega"}; | 
| 414 |  |  | } elsif ("$spec" eq "CIE-u") { | 
| 415 | greg | 2.54 | $cmd = q{rcalc -e "Ri=$1;Gi=$2;Bi=$3" } . | 
| 416 | greg | 2.53 | qq{-e "$CIEuv" } . | 
| 417 | greg | 2.54 | q{-e "$1=uprime"}; | 
| 418 | greg | 2.53 | } elsif ("$spec" eq "CIE-v") { | 
| 419 | greg | 2.54 | $cmd = q{rcalc -e "Ri=$1;Gi=$2;Bi=$3" } . | 
| 420 | greg | 2.53 | qq{-e "$CIEuv" } . | 
| 421 | greg | 2.54 | q{-e "$1=vprime"}; | 
| 422 | greg | 2.49 | } | 
| 423 | greg | 2.47 | } else { | 
| 424 | greg | 2.49 | if ("$spec" eq "Visible") { | 
| 425 |  |  | $cmd = "rcalc -if3 -e 'Omega:PI/($ns*$ns)' " . | 
| 426 | greg | 2.53 | q{-e 'Ri=$1;Gi=$2;Bi=$3' } . | 
| 427 |  |  | "-e '$CIEuv' " . | 
| 428 |  |  | q{-e '$1=Yi/Omega'}; | 
| 429 |  |  | } elsif ("$spec" eq "CIE-u") { | 
| 430 | greg | 2.54 | $cmd = q{rcalc -if3 -e 'Ri=$1;Gi=$2;Bi=$3' } . | 
| 431 | greg | 2.53 | "-e '$CIEuv' " . | 
| 432 | greg | 2.54 | q{-e '$1=uprime'}; | 
| 433 | greg | 2.53 | } elsif ("$spec" eq "CIE-v") { | 
| 434 | greg | 2.54 | $cmd = q{rcalc -if3 -e 'Ri=$1;Gi=$2;Bi=$3' } . | 
| 435 | greg | 2.53 | "-e '$CIEuv' " . | 
| 436 | greg | 2.54 | q{-e '$1=vprime'}; | 
| 437 | greg | 2.49 | } | 
| 438 | greg | 2.47 | } | 
| 439 | greg | 2.49 | if ($pctcull >= 0) { | 
| 440 | greg | 2.88 | my $avg = ( $dorecip && ( $tensortree == 3 || "$typ" =~ /^r[fb]/ ) ) ? " -a" : ""; | 
| 441 | greg | 2.53 | my $pcull = ("$spec" eq "Visible") ? $pctcull : | 
| 442 | greg | 2.55 | (100 - (100-$pctcull)*.25) ; | 
| 443 | greg | 2.49 | if ($windoz) { | 
| 444 |  |  | $cmd = "rcollate -ho -oc 1 $src | " . | 
| 445 |  |  | $cmd . | 
| 446 | greg | 2.53 | " | rttree_reduce$avg -h -fa -t $pcull -r $tensortree -g $ttlog2"; | 
| 447 | greg | 2.49 | } else { | 
| 448 |  |  | $cmd .= " -of $src " . | 
| 449 | greg | 2.53 | "| rttree_reduce$avg -h -ff -t $pcull -r $tensortree -g $ttlog2"; | 
| 450 | greg | 2.49 | } | 
| 451 | greg | 2.70 | run_check "$cmd > $dest"; | 
| 452 | greg | 2.47 | } else { | 
| 453 | greg | 2.49 | if ($windoz) { | 
| 454 |  |  | $cmd = "rcollate -ho -oc 1 $src | " . $cmd ; | 
| 455 |  |  | } else { | 
| 456 |  |  | $cmd .= " $src"; | 
| 457 |  |  | } | 
| 458 | greg | 2.70 | if ( active_phase() ) { | 
| 459 |  |  | open(DATOUT, "> $dest"); | 
| 460 |  |  | print DATOUT "{\n"; | 
| 461 |  |  | close DATOUT; | 
| 462 |  |  | # print STDERR "Running: $cmd\n"; | 
| 463 |  |  | system "$cmd >> $dest"; | 
| 464 |  |  | die "Failure running rcalc" if ( $? ); | 
| 465 |  |  | open(DATOUT, ">> $dest"); | 
| 466 |  |  | for (my $i = ($tensortree==3)*$ns*$ns*$ns/2; $i-- > 0; ) { | 
| 467 |  |  | print DATOUT "0\n"; | 
| 468 |  |  | } | 
| 469 |  |  | print DATOUT "}\n"; | 
| 470 |  |  | close DATOUT; | 
| 471 | greg | 2.49 | } | 
| 472 | greg | 2.47 | } | 
| 473 | greg | 2.49 | if ( "$spec" ne "$curspec" ) { | 
| 474 | greg | 2.51 | $wrapper .= " -s $spec"; | 
| 475 | greg | 2.49 | $curspec = $spec; | 
| 476 | greg | 2.34 | } | 
| 477 | greg | 2.49 | $wrapper .= " -$typ $dest"; | 
| 478 |  |  | }       # end of ttree_comp() | 
| 479 | greg | 2.15 |  | 
| 480 |  |  | #------------- End of do_tree_bsdf() & subroutines -------------# | 
| 481 |  |  |  | 
| 482 |  |  | #+++++++++++++++ Klems matrix BSDF generation +++++++++++++++# | 
| 483 |  |  | sub do_matrix_bsdf { | 
| 484 | greg | 2.49 |  | 
| 485 |  |  | # Run rfluxmtx processes to compute each side | 
| 486 |  |  | do_matrix_dir(0) if ( $doback ); | 
| 487 |  |  | do_matrix_dir(1) if ( $doforw ); | 
| 488 |  |  |  | 
| 489 |  |  | }       # end of sub do_matrix_bsdf() | 
| 490 |  |  |  | 
| 491 |  |  | # Call rfluxmtx and process tensor tree BSDF for the given direction | 
| 492 |  |  | sub do_matrix_dir { | 
| 493 |  |  | my $forw = shift; | 
| 494 | greg | 2.70 | my $dop = do_phase(); | 
| 495 |  |  | my $r = ($dop < 0) ? " -r" : ""; | 
| 496 | greg | 2.49 | my $sender = ($bsender,$fsender)[$forw]; | 
| 497 | greg | 2.71 | my $cmd = "$rfluxmtx$r -fd $sender $receivers -i $octree"; | 
| 498 |  |  | if ( $dop ) { | 
| 499 |  |  | # print STDERR "Running: $cmd\n"; | 
| 500 |  |  | system $cmd; | 
| 501 |  |  | die "Failure running rfluxmtx" if ( $? ); | 
| 502 |  |  | } | 
| 503 | greg | 2.49 | matrix_out($forw); | 
| 504 |  |  | }       # end of do_matrix_dir() | 
| 505 |  |  |  | 
| 506 |  |  | sub matrix_out { | 
| 507 |  |  | my $forw = shift; | 
| 508 |  |  | my ($refldat,$transdat); | 
| 509 |  |  | if ( $forw ) { | 
| 510 |  |  | $transdat = $facedat; | 
| 511 |  |  | $refldat = $behinddat; | 
| 512 | greg | 2.47 | } else { | 
| 513 | greg | 2.49 | $transdat = $behinddat; | 
| 514 |  |  | $refldat = $facedat; | 
| 515 | greg | 2.47 | } | 
| 516 | greg | 2.49 | # Output transmission | 
| 517 |  |  | my $ttyp = ("tb","tf")[$forw]; | 
| 518 |  |  | matrix_comp($ttyp, "Visible", $transdat, ($tb,$tf)[$forw]); | 
| 519 |  |  | if ( $docolor ) { | 
| 520 |  |  | matrix_comp($ttyp, "CIE-X", $transdat, ($tbx,$tfx)[$forw]); | 
| 521 |  |  | matrix_comp($ttyp, "CIE-Z", $transdat, ($tbz,$tfz)[$forw]); | 
| 522 |  |  | } | 
| 523 |  |  | # Output reflection | 
| 524 |  |  | my $rtyp = ("rb","rf")[$forw]; | 
| 525 |  |  | matrix_comp($rtyp, "Visible", $refldat, ($rb,$rf)[$forw]); | 
| 526 |  |  | if ( $docolor ) { | 
| 527 |  |  | matrix_comp($rtyp, "CIE-X", $refldat, ($rbx,$rfx)[$forw]); | 
| 528 |  |  | matrix_comp($rtyp, "CIE-Z", $refldat, ($rbz,$rfz)[$forw]); | 
| 529 |  |  | } | 
| 530 |  |  | }       # end of matrix_out() | 
| 531 |  |  |  | 
| 532 |  |  | # Transpose matrix component data and save to file | 
| 533 |  |  | sub matrix_comp { | 
| 534 |  |  | my $typ = shift; | 
| 535 |  |  | my $spec = shift; | 
| 536 |  |  | my $src = shift; | 
| 537 |  |  | my $dest = shift; | 
| 538 |  |  | my $cmd = "rmtxop -fa -t"; | 
| 539 |  |  | if ("$spec" eq "Visible") { | 
| 540 | greg | 2.55 | $cmd .= " -c 0.2651 0.6701 0.0648"; | 
| 541 | greg | 2.49 | } elsif ("$spec" eq "CIE-X") { | 
| 542 | greg | 2.55 | $cmd .= " -c 0.5141 0.3239 0.1620"; | 
| 543 | greg | 2.49 | } elsif ("$spec" eq "CIE-Z") { | 
| 544 | greg | 2.55 | $cmd .= " -c 0.0241 0.1229 0.8530"; | 
| 545 | greg | 2.49 | } | 
| 546 | greg | 2.72 | $cmd .= " $src | getinfo -"; | 
| 547 | greg | 2.70 | run_check "$cmd > $dest"; | 
| 548 | greg | 2.49 | if ( "$spec" ne "$curspec" ) { | 
| 549 | greg | 2.51 | $wrapper .= " -s $spec"; | 
| 550 | greg | 2.49 | $curspec = $spec; | 
| 551 | greg | 2.47 | } | 
| 552 | greg | 2.49 | $wrapper .= " -$typ $dest"; | 
| 553 |  |  | }       # end of matrix_comp() | 
| 554 | greg | 2.47 |  | 
| 555 | greg | 2.49 | #------------- End of do_matrix_bsdf() & subroutines --------------# |