ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/util/genBSDF.pl
Revision: 2.85
Committed: Tue May 5 22:35:58 2020 UTC (3 years, 11 months ago) by greg
Content type: text/plain
Branch: MAIN
Changes since 2.84: +3 -3 lines
Log Message:
Removed unnecessary (and incorrect) isotropic sampling constraint

File Contents

# User Rev Content
1 greg 2.1 #!/usr/bin/perl -w
2 greg 2.85 # RCSid $Id: genBSDF.pl,v 2.84 2020/01/17 17:33:48 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     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};
217     } else {
218     system "rad2mgf $radscn >> $mgfscn";
219     }
220     open(MGFSCN, ">> $mgfscn");
221     print MGFSCN "xf\n";
222     close MGFSCN;
223     }
224     # Create receiver & sender surfaces (rectangular)
225     open(RADSCN, "> $receivers");
226     print RADSCN '#@rfluxmtx ' . ($tensortree ? "h=-sc$ns\n" : "h=-kf\n");
227     print RADSCN '#@rfluxmtx ' . "u=-Y o=$facedat\n\n";
228     print RADSCN "void glow receiver_face\n0\n0\n4 1 1 1 0\n\n";
229     print RADSCN "receiver_face source f_receiver\n0\n0\n4 0 0 1 180\n\n";
230     print RADSCN '#@rfluxmtx ' . ($tensortree ? "h=+sc$ns\n" : "h=+kf\n");
231     print RADSCN '#@rfluxmtx ' . "u=-Y o=$behinddat\n\n";
232     print RADSCN "void glow receiver_behind\n0\n0\n4 1 1 1 0\n\n";
233     print RADSCN "receiver_behind source b_receiver\n0\n0\n4 0 0 -1 180\n";
234 greg 2.49 close RADSCN;
235 greg 2.66 # Prepare sender surfaces
236     if ( $tensortree != 3 ) { # Isotropic tensor tree is exception
237     open (RADSCN, "> $fsender");
238     print RADSCN '#@rfluxmtx u=-Y ' . ($tensortree ? "h=-sc$ns\n\n" : "h=-kf\n\n");
239     print RADSCN "void polygon fwd_sender\n0\n0\n12\n";
240     printf RADSCN "\t%e\t%e\t%e\n", $dim[0], $dim[2], $dim[4]-$FEPS;
241     printf RADSCN "\t%e\t%e\t%e\n", $dim[0], $dim[3], $dim[4]-$FEPS;
242     printf RADSCN "\t%e\t%e\t%e\n", $dim[1], $dim[3], $dim[4]-$FEPS;
243     printf RADSCN "\t%e\t%e\t%e\n", $dim[1], $dim[2], $dim[4]-$FEPS;
244     close RADSCN;
245     open (RADSCN, "> $bsender");
246     print RADSCN '#@rfluxmtx u=-Y ' . ($tensortree ? "h=+sc$ns\n\n" : "h=+kf\n\n");
247     print RADSCN "void polygon bwd_sender\n0\n0\n12\n";
248     printf RADSCN "\t%e\t%e\t%e\n", $dim[0], $dim[2], $dim[5]+$FEPS;
249     printf RADSCN "\t%e\t%e\t%e\n", $dim[1], $dim[2], $dim[5]+$FEPS;
250     printf RADSCN "\t%e\t%e\t%e\n", $dim[1], $dim[3], $dim[5]+$FEPS;
251     printf RADSCN "\t%e\t%e\t%e\n", $dim[0], $dim[3], $dim[5]+$FEPS;
252     close RADSCN;
253     }
254     print STDERR "Recover using: $0 -recover $td\n";
255     }
256     $curphase = 0;
257 greg 2.15 # Create data segments (all the work happens here)
258     if ( $tensortree ) {
259     do_tree_bsdf();
260     } else {
261     do_matrix_bsdf();
262     }
263 greg 2.49 # Output XML
264 greg 2.67 # print STDERR "Running: $wrapper\n";
265 greg 2.62 system "$wrapper -C \"Created by: genBSDF @savedARGV\"";
266 greg 2.49 die "Could not wrap BSDF data\n" if ( $? );
267 greg 2.15 # Clean up temporary files and exit
268 greg 2.53 exec $rmtmp;
269 greg 2.16
270 greg 2.67 #============== End of main program segment ==============#
271    
272     # Function to determine if next phase should be skipped or recovered
273     sub do_phase {
274     $curphase++;
275 greg 2.70 if ( defined $recovery ) {
276     if ( $recovery > $curphase ) { return 0; }
277     if ( $recovery == $curphase ) { return -1; }
278 greg 2.67 }
279 greg 2.73 open(MYPH, ">> $td/phase.txt");
280 greg 2.67 print MYPH "$curphase\n";
281 greg 2.73 close MYPH;
282 greg 2.67 return 1;
283     }
284 greg 2.15
285 greg 2.70 # Check if we are in active phase (not skipping parts)
286     sub active_phase {
287     if ( defined $recovery ) {
288     if ( $recovery > $curphase ) { return 0; }
289     if ( $recovery == $curphase ) { return -1; }
290     }
291     return 1;
292     }
293    
294     # Function to run program and check output if in active phase
295     sub run_check {
296     if ( !active_phase() ) { return; }
297     my $cmd = shift;
298     # print STDERR "Running: $cmd\n";
299     system $cmd;
300     die "Failure running: $cmd\n" if ( $? );
301     }
302    
303 greg 2.15 #++++++++++++++ Tensor tree BSDF generation ++++++++++++++#
304     sub do_tree_bsdf {
305 greg 2.49
306     # Run rfluxmtx processes to compute each side
307     do_ttree_dir(0) if ( $doback );
308     do_ttree_dir(1) if ( $doforw );
309 greg 2.36
310 greg 2.15 } # end of sub do_tree_bsdf()
311    
312 greg 2.49 # Call rfluxmtx and process tensor tree BSDF for the given direction
313     sub do_ttree_dir {
314 greg 2.15 my $forw = shift;
315 greg 2.70 my $dop = do_phase();
316     my $r = ($dop < 0) ? " -r" : "";
317 greg 2.47 my $cmd;
318 greg 2.15 if ( $tensortree == 3 ) {
319     # Isotropic BSDF
320 greg 2.36 my $ns2 = $ns / 2;
321 greg 2.47 if ($windoz) {
322     $cmd = "cnt $ns2 $ny $nx " .
323     qq{| rcalc -e "r1=rand(.8681*recno-.673892)" } .
324     qq{-e "r2=rand(-5.37138*recno+67.1737811)" } .
325     qq{-e "r3=rand(+3.17603772*recno+83.766771)" } .
326 greg 2.82 qq{-e "r4=rand(-1.5839226*recno-59.82712)" } .
327 greg 2.84 qq{-e "odds(n):if(.5*n-floor(.5*n)-.25,-1,1)"} .
328     qq{-e "Dx=1-(\$1+r1)/$ns2" } .
329 greg 2.85 qq{-e "Dy=odds(\$1)/$ns*r2" } .
330 greg 2.82 qq{-e "Dz=sqrt(1-Dx*Dx-Dy*Dy)" } .
331 greg 2.47 qq{-e "xp=(\$3+r2)*(($dim[1]-$dim[0])/$nx)+$dim[0]" } .
332     qq{-e "yp=(\$2+r3)*(($dim[3]-$dim[2])/$ny)+$dim[2]" } .
333     qq{-e "zp=$dim[5-$forw]" -e "myDz=Dz*($forw*2-1)" } .
334     qq{-e "\$1=xp-Dx;\$2=yp-Dy;\$3=zp-myDz" } .
335     qq{-e "\$4=Dx;\$5=Dy;\$6=myDz" } .
336 greg 2.66 "| $rfluxmtx$r -fa -y $ns2 - $receivers -i $octree";
337 greg 2.47 } else {
338     $cmd = "cnt $ns2 $ny $nx " .
339 greg 2.75 qq{| rcalc -e 'r1=rand(.8681*recno-.673892)' } .
340     qq{-e 'r2=rand(-5.37138*recno+67.1737811)' } .
341     qq{-e 'r3=rand(+3.17603772*recno+83.766771)' } .
342 greg 2.82 qq{-e 'r4=rand(-1.5839226*recno-59.82712)' } .
343 greg 2.84 qq{-e 'odds(n):if(.5*n-floor(.5*n)-.25,-1,1)' } .
344     qq{-e 'Dx=1-(\$1+r1)/$ns2' } .
345 greg 2.85 qq{-e 'Dy=odds(\$1)/$ns*r2' } .
346 greg 2.82 qq{-e 'Dz=sqrt(1-Dx*Dx-Dy*Dy)' } .
347     qq{-e 'xp=(\$3+r3)*(($dim[1]-$dim[0])/$nx)+$dim[0]' } .
348     qq{-e 'yp=(\$2+r4)*(($dim[3]-$dim[2])/$ny)+$dim[2]' } .
349 greg 2.75 qq{-e 'zp=$dim[5-$forw]' -e 'myDz=Dz*($forw*2-1)' } .
350 greg 2.47 qq{-e '\$1=xp-Dx;\$2=yp-Dy;\$3=zp-myDz' } .
351     qq{-e '\$4=Dx;\$5=Dy;\$6=myDz' -of } .
352 greg 2.66 "| $rfluxmtx$r -h -ff -y $ns2 - $receivers -i $octree";
353 greg 2.47 }
354 greg 2.15 } else {
355     # Anisotropic BSDF
356 greg 2.49 my $sender = ($bsender,$fsender)[$forw];
357 greg 2.47 if ($windoz) {
358 greg 2.66 $cmd = "$rfluxmtx$r -fa $sender $receivers -i $octree";
359 greg 2.47 } else {
360 greg 2.66 $cmd = "$rfluxmtx$r -h -ff $sender $receivers -i $octree";
361 greg 2.47 }
362 greg 2.15 }
363 greg 2.71 if ( $dop ) {
364     # print STDERR "Running: $cmd\n";
365     system $cmd;
366     die "Failure running rfluxmtx" if ( $? );
367     }
368 greg 2.36 ttree_out($forw);
369 greg 2.49 } # end of do_ttree_dir()
370 greg 2.15
371 greg 2.49 # Simplify and store tensor tree results
372 greg 2.15 sub ttree_out {
373     my $forw = shift;
374 greg 2.49 my ($refldat,$transdat);
375     if ( $forw ) {
376     $transdat = $facedat;
377     $refldat = $behinddat;
378 greg 2.47 } else {
379 greg 2.49 $transdat = $behinddat;
380     $refldat = $facedat;
381 greg 2.47 }
382 greg 2.49 # Only output one transmitted anisotropic distribution, preferring backwards
383     if ( !$forw || !$doback || $tensortree==3 ) {
384     my $ttyp = ("tb","tf")[$forw];
385     ttree_comp($ttyp, "Visible", $transdat, ($tb,$tf)[$forw]);
386     if ( $docolor ) {
387 greg 2.53 ttree_comp($ttyp, "CIE-u", $transdat, ($tbx,$tfx)[$forw]);
388     ttree_comp($ttyp, "CIE-v", $transdat, ($tbz,$tfz)[$forw]);
389 greg 2.49 }
390 greg 2.47 }
391 greg 2.49 # Output reflection
392     my $rtyp = ("rb","rf")[$forw];
393     ttree_comp($rtyp, "Visible", $refldat, ($rb,$rf)[$forw]);
394     if ( $docolor ) {
395 greg 2.53 ttree_comp($rtyp, "CIE-u", $refldat, ($rbx,$rfx)[$forw]);
396     ttree_comp($rtyp, "CIE-v", $refldat, ($rbz,$rfz)[$forw]);
397 greg 2.34 }
398 greg 2.49 } # end of ttree_out()
399    
400     # Call rttree_reduce on the given component
401     sub ttree_comp {
402     my $typ = shift;
403     my $spec = shift;
404     my $src = shift;
405     my $dest = shift;
406     my $cmd;
407 greg 2.47 if ($windoz) {
408 greg 2.49 if ("$spec" eq "Visible") {
409     $cmd = qq{rcalc -e "Omega:PI/($ns*$ns)" } .
410 greg 2.53 q{-e "Ri=$1;Gi=$2;Bi=$3" } .
411     qq{-e "$CIEuv" } .
412     q{-e "$1=Yi/Omega"};
413     } elsif ("$spec" eq "CIE-u") {
414 greg 2.54 $cmd = q{rcalc -e "Ri=$1;Gi=$2;Bi=$3" } .
415 greg 2.53 qq{-e "$CIEuv" } .
416 greg 2.54 q{-e "$1=uprime"};
417 greg 2.53 } elsif ("$spec" eq "CIE-v") {
418 greg 2.54 $cmd = q{rcalc -e "Ri=$1;Gi=$2;Bi=$3" } .
419 greg 2.53 qq{-e "$CIEuv" } .
420 greg 2.54 q{-e "$1=vprime"};
421 greg 2.49 }
422 greg 2.47 } else {
423 greg 2.49 if ("$spec" eq "Visible") {
424     $cmd = "rcalc -if3 -e 'Omega:PI/($ns*$ns)' " .
425 greg 2.53 q{-e 'Ri=$1;Gi=$2;Bi=$3' } .
426     "-e '$CIEuv' " .
427     q{-e '$1=Yi/Omega'};
428     } elsif ("$spec" eq "CIE-u") {
429 greg 2.54 $cmd = q{rcalc -if3 -e 'Ri=$1;Gi=$2;Bi=$3' } .
430 greg 2.53 "-e '$CIEuv' " .
431 greg 2.54 q{-e '$1=uprime'};
432 greg 2.53 } elsif ("$spec" eq "CIE-v") {
433 greg 2.54 $cmd = q{rcalc -if3 -e 'Ri=$1;Gi=$2;Bi=$3' } .
434 greg 2.53 "-e '$CIEuv' " .
435 greg 2.54 q{-e '$1=vprime'};
436 greg 2.49 }
437 greg 2.47 }
438 greg 2.49 if ($pctcull >= 0) {
439 greg 2.76 my $avg = ( $dorecip && "$typ" =~ /^r[fb]/ ) ? " -a" : "";
440 greg 2.53 my $pcull = ("$spec" eq "Visible") ? $pctcull :
441 greg 2.55 (100 - (100-$pctcull)*.25) ;
442 greg 2.49 if ($windoz) {
443     $cmd = "rcollate -ho -oc 1 $src | " .
444     $cmd .
445 greg 2.53 " | rttree_reduce$avg -h -fa -t $pcull -r $tensortree -g $ttlog2";
446 greg 2.49 } else {
447     $cmd .= " -of $src " .
448 greg 2.53 "| rttree_reduce$avg -h -ff -t $pcull -r $tensortree -g $ttlog2";
449 greg 2.49 }
450 greg 2.70 run_check "$cmd > $dest";
451 greg 2.47 } else {
452 greg 2.49 if ($windoz) {
453     $cmd = "rcollate -ho -oc 1 $src | " . $cmd ;
454     } else {
455     $cmd .= " $src";
456     }
457 greg 2.70 if ( active_phase() ) {
458     open(DATOUT, "> $dest");
459     print DATOUT "{\n";
460     close DATOUT;
461     # print STDERR "Running: $cmd\n";
462     system "$cmd >> $dest";
463     die "Failure running rcalc" if ( $? );
464     open(DATOUT, ">> $dest");
465     for (my $i = ($tensortree==3)*$ns*$ns*$ns/2; $i-- > 0; ) {
466     print DATOUT "0\n";
467     }
468     print DATOUT "}\n";
469     close DATOUT;
470 greg 2.49 }
471 greg 2.47 }
472 greg 2.49 if ( "$spec" ne "$curspec" ) {
473 greg 2.51 $wrapper .= " -s $spec";
474 greg 2.49 $curspec = $spec;
475 greg 2.34 }
476 greg 2.49 $wrapper .= " -$typ $dest";
477     } # end of ttree_comp()
478 greg 2.15
479     #------------- End of do_tree_bsdf() & subroutines -------------#
480    
481     #+++++++++++++++ Klems matrix BSDF generation +++++++++++++++#
482     sub do_matrix_bsdf {
483 greg 2.49
484     # Run rfluxmtx processes to compute each side
485     do_matrix_dir(0) if ( $doback );
486     do_matrix_dir(1) if ( $doforw );
487    
488     } # end of sub do_matrix_bsdf()
489    
490     # Call rfluxmtx and process tensor tree BSDF for the given direction
491     sub do_matrix_dir {
492     my $forw = shift;
493 greg 2.70 my $dop = do_phase();
494     my $r = ($dop < 0) ? " -r" : "";
495 greg 2.49 my $sender = ($bsender,$fsender)[$forw];
496 greg 2.71 my $cmd = "$rfluxmtx$r -fd $sender $receivers -i $octree";
497     if ( $dop ) {
498     # print STDERR "Running: $cmd\n";
499     system $cmd;
500     die "Failure running rfluxmtx" if ( $? );
501     }
502 greg 2.49 matrix_out($forw);
503     } # end of do_matrix_dir()
504    
505     sub matrix_out {
506     my $forw = shift;
507     my ($refldat,$transdat);
508     if ( $forw ) {
509     $transdat = $facedat;
510     $refldat = $behinddat;
511 greg 2.47 } else {
512 greg 2.49 $transdat = $behinddat;
513     $refldat = $facedat;
514 greg 2.47 }
515 greg 2.49 # Output transmission
516     my $ttyp = ("tb","tf")[$forw];
517     matrix_comp($ttyp, "Visible", $transdat, ($tb,$tf)[$forw]);
518     if ( $docolor ) {
519     matrix_comp($ttyp, "CIE-X", $transdat, ($tbx,$tfx)[$forw]);
520     matrix_comp($ttyp, "CIE-Z", $transdat, ($tbz,$tfz)[$forw]);
521     }
522     # Output reflection
523     my $rtyp = ("rb","rf")[$forw];
524     matrix_comp($rtyp, "Visible", $refldat, ($rb,$rf)[$forw]);
525     if ( $docolor ) {
526     matrix_comp($rtyp, "CIE-X", $refldat, ($rbx,$rfx)[$forw]);
527     matrix_comp($rtyp, "CIE-Z", $refldat, ($rbz,$rfz)[$forw]);
528     }
529     } # end of matrix_out()
530    
531     # Transpose matrix component data and save to file
532     sub matrix_comp {
533     my $typ = shift;
534     my $spec = shift;
535     my $src = shift;
536     my $dest = shift;
537     my $cmd = "rmtxop -fa -t";
538     if ("$spec" eq "Visible") {
539 greg 2.55 $cmd .= " -c 0.2651 0.6701 0.0648";
540 greg 2.49 } elsif ("$spec" eq "CIE-X") {
541 greg 2.55 $cmd .= " -c 0.5141 0.3239 0.1620";
542 greg 2.49 } elsif ("$spec" eq "CIE-Z") {
543 greg 2.55 $cmd .= " -c 0.0241 0.1229 0.8530";
544 greg 2.49 }
545 greg 2.72 $cmd .= " $src | getinfo -";
546 greg 2.70 run_check "$cmd > $dest";
547 greg 2.49 if ( "$spec" ne "$curspec" ) {
548 greg 2.51 $wrapper .= " -s $spec";
549 greg 2.49 $curspec = $spec;
550 greg 2.47 }
551 greg 2.49 $wrapper .= " -$typ $dest";
552     } # end of matrix_comp()
553 greg 2.47
554 greg 2.49 #------------- End of do_matrix_bsdf() & subroutines --------------#