ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/util/genBSDF.pl
Revision: 2.91
Committed: Fri Jan 5 18:04:28 2024 UTC (3 months, 3 weeks ago) by greg
Content type: text/plain
Branch: MAIN
CVS Tags: HEAD
Changes since 2.90: +2 -2 lines
Log Message:
fix(genBSDF): Enforce the fact that spectral sampling is currently unsupported

File Contents

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