ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/util/genBSDF.pl
Revision: 2.89
Committed: Wed Mar 2 21:15:18 2022 UTC (2 years ago) by greg
Content type: text/plain
Branch: MAIN
Changes since 2.88: +2 -1 lines
Log Message:
fix(genBSDF): Added error message for +mgf with +geom and no input file

File Contents

# Content
1 #!/usr/bin/perl -w
2 # RCSid $Id: genBSDF.pl,v 2.88 2022/01/24 17:49: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,$cph);
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 $cph = "$td\\phase.txt";
61 $rmtmp = "rd /S /Q $td";
62 } else {
63 $radscn = "$td/device.rad";
64 $mgfscn = "$td/device.mgf";
65 $octree = "$td/device.oct";
66 $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 $cph = "$td/phase.txt";
84 $rmtmp = "rm -rf $td";
85 }
86 my @savedARGV = @ARGV;
87 my $rfluxmtx = "rfluxmtx -ab 5 -ad 700 -lw 3e-6 -w-";
88 my $wrapper = "wrapBSDF";
89 my $tensortree = 0;
90 my $ttlog2 = 4;
91 my $dorecip = 1;
92 my $nsamp = 2000;
93 my $mgfin = 0;
94 my $geout = 1;
95 my $nproc = 1;
96 my $docolor = 0;
97 my $doforw = 0;
98 my $doback = 1;
99 my $pctcull = 90;
100 my $gunit = "meter";
101 my $curspec = "Visible";
102 my @dim;
103 # Get options
104 while ($#ARGV >= 0) {
105 if ("$ARGV[0]" =~ /^[-+]m/) {
106 $mgfin = ("$ARGV[0]" =~ /^\+/);
107 } elsif ("$ARGV[0]" eq "-r") {
108 $rfluxmtx .= " $ARGV[1]";
109 shift @ARGV;
110 } elsif ("$ARGV[0]" =~ /^[-+]g/) {
111 $geout = ("$ARGV[0]" =~ /^\+/);
112 $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 } elsif ("$ARGV[0]" =~ /^[-+]C/) {
118 $docolor = ("$ARGV[0]" =~ /^\+/);
119 } elsif ("$ARGV[0]" =~ /^[-+]a/) {
120 $dorecip = ("$ARGV[0]" =~ /^\+/);
121 } elsif ("$ARGV[0]" =~ /^[-+]f/) {
122 $doforw = ("$ARGV[0]" =~ /^\+/);
123 } elsif ("$ARGV[0]" =~ /^[-+]b/) {
124 $doback = ("$ARGV[0]" =~ /^\+/);
125 } elsif ("$ARGV[0]" eq "-t") {
126 # Use value < 0 for rttree_reduce bypass
127 $pctcull = $ARGV[1];
128 if ($pctcull >= 100) {
129 die "Illegal -t culling percentage, must be < 100\n";
130 }
131 shift @ARGV;
132 } elsif ("$ARGV[0]" =~ /^-t[34]$/) {
133 $tensortree = substr($ARGV[0], 2, 1);
134 $ttlog2 = $ARGV[1];
135 shift @ARGV;
136 } elsif ("$ARGV[0]" eq "-s") {
137 $wrapper .= " -f \"$ARGV[1]\"";
138 shift @ARGV;
139 } elsif ("$ARGV[0]" eq "-W") {
140 $wrapper .= " -W";
141 } 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 @dim = @ARGV[1..6];
150 shift @ARGV for (1..6);
151 } elsif ("$ARGV[0]" =~ /^[-+]./) {
152 userror();
153 } else {
154 last;
155 }
156 shift @ARGV;
157 }
158 # Check that we're actually being asked to do something
159 die "Must have at least one of +forward or +backward\n" if (!$doforw && !$doback);
160 $wrapper .= $tensortree ? " -a t$tensortree" : " -a kf -c";
161 $wrapper .= " -u $gunit";
162 if ( !defined $recovery ) {
163 # 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 # Get scene description
167 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 }
175 if ( $#dim != 5 ) {
176 @dim = split ' ', `getbbox -h -w $radscn`;
177 }
178 die "Device entirely inside room!\n" if ( $dim[4] >= 0 );
179 if ( $dim[5] > 1e-5 ) {
180 print STDERR "Warning: Device extends into room\n";
181 } elsif ( $dim[5]*$dim[5] > .01*($dim[1]-$dim[0])*($dim[3]-$dim[2]) ) {
182 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 ';h=' . ($dim[3] - $dim[2]) . '"';
187 $wrapper .= " -g $mgfscn" if ( $geout );
188 # 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 '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 my $FEPS = 1e-5;
196 my $ns = 2**$ttlog2;
197 my $nx = int(sqrt($nsamp*($dim[1]-$dim[0])/($dim[3]-$dim[2])) + 1);
198 my $ny = int($nsamp/$nx + 1);
199 $nsamp = $nx * $ny;
200 $rfluxmtx .= " -n $nproc -c $nsamp";
201 if ( !defined $recovery ) {
202 open(MYAVH, "> $td/savedARGV.txt");
203 foreach (@savedARGV) {
204 print MYAVH "$_\n";
205 }
206 close MYAVH;
207 # Generate octree
208 system "oconv -w $radscn > $octree";
209 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 die "+mgf requires input file with +geom\n" if ($#ARGV < 0);
217 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 close RADSCN;
236 # 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 # Create data segments (all the work happens here)
259 if ( $tensortree ) {
260 do_tree_bsdf();
261 } else {
262 do_matrix_bsdf();
263 }
264 # Output XML
265 # print STDERR "Running: $wrapper\n";
266 system "$wrapper -C \"Created by: genBSDF @savedARGV\"";
267 die "Could not wrap BSDF data\n" if ( $? );
268 # Clean up temporary files and exit
269 exec $rmtmp;
270
271 #============== 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 if ( defined $recovery ) {
277 if ( $recovery > $curphase ) { return 0; }
278 if ( $recovery == $curphase ) { return -1; }
279 }
280 open(MYPH, ">> $td/phase.txt");
281 print MYPH "$curphase\n";
282 close MYPH;
283 return 1;
284 }
285
286 # 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 #++++++++++++++ Tensor tree BSDF generation ++++++++++++++#
305 sub do_tree_bsdf {
306
307 # Run rfluxmtx processes to compute each side
308 do_ttree_dir(0) if ( $doback );
309 do_ttree_dir(1) if ( $doforw );
310
311 } # end of sub do_tree_bsdf()
312
313 # Call rfluxmtx and process tensor tree BSDF for the given direction
314 sub do_ttree_dir {
315 my $forw = shift;
316 my $dop = do_phase();
317 my $r = ($dop < 0) ? " -r" : "";
318 my $cmd;
319 if ( $tensortree == 3 ) {
320 # Isotropic BSDF
321 my $ns2 = $ns / 2;
322 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 qq{-e "r4=rand(-1.5839226*recno-59.82712)" } .
328 qq{-e "odds(n):if(.5*n-floor(.5*n)-.25,-1,1)"} .
329 qq{-e "Dx=1-(\$1+r1)/$ns2" } .
330 qq{-e "Dy=min(1/$ns,sqrt(1-Dx*Dx))*odds(\$1)*r2" } .
331 qq{-e "Dz=sqrt(1-Dx*Dx-Dy*Dy)" } .
332 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 "| $rfluxmtx$r -fa -y $ns2 - $receivers -i $octree";
338 } else {
339 $cmd = "cnt $ns2 $ny $nx " .
340 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 qq{-e 'r4=rand(-1.5839226*recno-59.82712)' } .
344 qq{-e 'odds(n):if(.5*n-floor(.5*n)-.25,-1,1)' } .
345 qq{-e 'Dx=1-(\$1+r1)/$ns2' } .
346 qq{-e 'Dy=min(1/$ns,sqrt(1-Dx*Dx))*odds(\$1)*r2' } .
347 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 qq{-e 'zp=$dim[5-$forw]' -e 'myDz=Dz*($forw*2-1)' } .
351 qq{-e '\$1=xp-Dx;\$2=yp-Dy;\$3=zp-myDz' } .
352 qq{-e '\$4=Dx;\$5=Dy;\$6=myDz' -of } .
353 "| $rfluxmtx$r -h -ff -y $ns2 - $receivers -i $octree";
354 }
355 } else {
356 # Anisotropic BSDF
357 my $sender = ($bsender,$fsender)[$forw];
358 if ($windoz) {
359 $cmd = "$rfluxmtx$r -fa $sender $receivers -i $octree";
360 } else {
361 $cmd = "$rfluxmtx$r -h -ff $sender $receivers -i $octree";
362 }
363 }
364 if ( $dop ) {
365 # print STDERR "Running: $cmd\n";
366 system $cmd;
367 die "Failure running rfluxmtx" if ( $? );
368 }
369 ttree_out($forw);
370 } # end of do_ttree_dir()
371
372 # Simplify and store tensor tree results
373 sub ttree_out {
374 my $forw = shift;
375 my ($refldat,$transdat);
376 if ( $forw ) {
377 $transdat = $facedat;
378 $refldat = $behinddat;
379 } else {
380 $transdat = $behinddat;
381 $refldat = $facedat;
382 }
383 # 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 ttree_comp($ttyp, "CIE-u", $transdat, ($tbx,$tfx)[$forw]);
389 ttree_comp($ttyp, "CIE-v", $transdat, ($tbz,$tfz)[$forw]);
390 }
391 }
392 # Output reflection
393 my $rtyp = ("rb","rf")[$forw];
394 ttree_comp($rtyp, "Visible", $refldat, ($rb,$rf)[$forw]);
395 if ( $docolor ) {
396 ttree_comp($rtyp, "CIE-u", $refldat, ($rbx,$rfx)[$forw]);
397 ttree_comp($rtyp, "CIE-v", $refldat, ($rbz,$rfz)[$forw]);
398 }
399 } # 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 if ($windoz) {
409 if ("$spec" eq "Visible") {
410 $cmd = qq{rcalc -e "Omega:PI/($ns*$ns)" } .
411 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 $cmd = q{rcalc -e "Ri=$1;Gi=$2;Bi=$3" } .
416 qq{-e "$CIEuv" } .
417 q{-e "$1=uprime"};
418 } elsif ("$spec" eq "CIE-v") {
419 $cmd = q{rcalc -e "Ri=$1;Gi=$2;Bi=$3" } .
420 qq{-e "$CIEuv" } .
421 q{-e "$1=vprime"};
422 }
423 } else {
424 if ("$spec" eq "Visible") {
425 $cmd = "rcalc -if3 -e 'Omega:PI/($ns*$ns)' " .
426 q{-e 'Ri=$1;Gi=$2;Bi=$3' } .
427 "-e '$CIEuv' " .
428 q{-e '$1=Yi/Omega'};
429 } elsif ("$spec" eq "CIE-u") {
430 $cmd = q{rcalc -if3 -e 'Ri=$1;Gi=$2;Bi=$3' } .
431 "-e '$CIEuv' " .
432 q{-e '$1=uprime'};
433 } elsif ("$spec" eq "CIE-v") {
434 $cmd = q{rcalc -if3 -e 'Ri=$1;Gi=$2;Bi=$3' } .
435 "-e '$CIEuv' " .
436 q{-e '$1=vprime'};
437 }
438 }
439 if ($pctcull >= 0) {
440 my $avg = ( $dorecip && ( $tensortree == 3 || "$typ" =~ /^r[fb]/ ) ) ? " -a" : "";
441 my $pcull = ("$spec" eq "Visible") ? $pctcull :
442 (100 - (100-$pctcull)*.25) ;
443 if ($windoz) {
444 $cmd = "rcollate -ho -oc 1 $src | " .
445 $cmd .
446 " | rttree_reduce$avg -h -fa -t $pcull -r $tensortree -g $ttlog2";
447 } else {
448 $cmd .= " -of $src " .
449 "| rttree_reduce$avg -h -ff -t $pcull -r $tensortree -g $ttlog2";
450 }
451 run_check "$cmd > $dest";
452 } else {
453 if ($windoz) {
454 $cmd = "rcollate -ho -oc 1 $src | " . $cmd ;
455 } else {
456 $cmd .= " $src";
457 }
458 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 }
472 }
473 if ( "$spec" ne "$curspec" ) {
474 $wrapper .= " -s $spec";
475 $curspec = $spec;
476 }
477 $wrapper .= " -$typ $dest";
478 } # end of ttree_comp()
479
480 #------------- End of do_tree_bsdf() & subroutines -------------#
481
482 #+++++++++++++++ Klems matrix BSDF generation +++++++++++++++#
483 sub do_matrix_bsdf {
484
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 my $dop = do_phase();
495 my $r = ($dop < 0) ? " -r" : "";
496 my $sender = ($bsender,$fsender)[$forw];
497 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 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 } else {
513 $transdat = $behinddat;
514 $refldat = $facedat;
515 }
516 # 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 $cmd .= " -c 0.2651 0.6701 0.0648";
541 } elsif ("$spec" eq "CIE-X") {
542 $cmd .= " -c 0.5141 0.3239 0.1620";
543 } elsif ("$spec" eq "CIE-Z") {
544 $cmd .= " -c 0.0241 0.1229 0.8530";
545 }
546 $cmd .= " $src | getinfo -";
547 run_check "$cmd > $dest";
548 if ( "$spec" ne "$curspec" ) {
549 $wrapper .= " -s $spec";
550 $curspec = $spec;
551 }
552 $wrapper .= " -$typ $dest";
553 } # end of matrix_comp()
554
555 #------------- End of do_matrix_bsdf() & subroutines --------------#