ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/util/genBSDF.pl
(Generate patch)

Comparing ray/src/util/genBSDF.pl (file contents):
Revision 2.3 by greg, Thu Sep 9 06:02:16 2010 UTC vs.
Revision 2.13 by greg, Sat Apr 16 00:39:07 2011 UTC

# Line 6 | Line 6
6   #       G. Ward
7   #
8   use strict;
9 + use File::Temp qw/ :mktemp  /;
10   sub userror {
11 <        print STDERR "Usage: genBSDF [-n Nproc][-c Nsamp][-dim xmin xmax ymin ymax zmin zmax][{+|-}mgf][{+|-}geom] [input ..]\n";
11 >        print STDERR "Usage: genBSDF [-n Nproc][-c Nsamp][-r \"ropts\"][-dim xmin xmax ymin ymax zmin zmax][{+|-}f][{+|-}b][{+|-}mgf][{+|-}geom] [input ..]\n";
12          exit 1;
13   }
14 < my $td = `mktemp -d /tmp/genBSDF.XXXXXX`;
14 > my $td = mkdtemp("/tmp/genBSDF.XXXXXX");
15   chomp $td;
16   my $nsamp = 1000;
17 + my $rtargs = "-w -ab 5 -ad 700 -lw 3e-6";
18   my $mgfin = 0;
19   my $geout = 1;
20   my $nproc = 1;
21 + my $doforw = 0;
22 + my $doback = 1;
23   my @dim;
24   # Get options
25   while ($#ARGV >= 0) {
26          if ("$ARGV[0]" =~ /^[-+]m/) {
27                  $mgfin = ("$ARGV[0]" =~ /^\+/);
28 +        } elsif ("$ARGV[0]" eq "-r") {
29 +                $rtargs = "$rtargs $ARGV[1]";
30 +                shift @ARGV;
31          } elsif ("$ARGV[0]" =~ /^[-+]g/) {
32                  $geout = ("$ARGV[0]" =~ /^\+/);
33 +        } elsif ("$ARGV[0]" =~ /^[-+]f/) {
34 +                $doforw = ("$ARGV[0]" =~ /^\+/);
35 +        } elsif ("$ARGV[0]" =~ /^[-+]b/) {
36 +                $doback = ("$ARGV[0]" =~ /^\+/);
37          } elsif ("$ARGV[0]" eq "-c") {
38                  $nsamp = $ARGV[1];
39                  shift @ARGV;
# Line 31 | Line 42 | while ($#ARGV >= 0) {
42                  shift @ARGV;
43          } elsif ("$ARGV[0]" =~ /^-d/) {
44                  userror() if ($#ARGV < 6);
45 <                @dim = "@ARGV[1..6]";
45 >                @dim = @ARGV[1..6];
46                  shift @ARGV for (1..6);
47          } elsif ("$ARGV[0]" =~ /^[-+]./) {
48                  userror();
# Line 40 | Line 51 | while ($#ARGV >= 0) {
51          }
52          shift @ARGV;
53   }
54 + # Check that we're actually being asked to do something
55 + die "Must have at least one of +forward or +backward" if (!$doforw && !$doback);
56   # Get scene description and dimensions
57   my $radscn = "$td/device.rad";
58   my $mgfscn = "$td/device.mgf";
# Line 49 | Line 62 | if ( $mgfin ) {
62          die "Could not load MGF input\n" if ( $? );
63          system "mgf2rad $mgfscn > $radscn";
64   } else {
65 <        system "cat @ARGV | xform -e > $radscn";
65 >        system "xform -e @ARGV > $radscn";
66          die "Could not load Radiance input\n" if ( $? );
67          system "rad2mgf $radscn > $mgfscn" if ( $geout );
68   }
69   if ($#dim != 5) {
70 <        @dim = split /\s+/, `getbbox -h $radscn`;
58 <        shift @dim;
70 >        @dim = split ' ', `getbbox -h $radscn`;
71   }
72   print STDERR "Warning: Device extends into room\n" if ($dim[5] > 1e-5);
73 < # Add receiver surface (rectangle)
74 < my $modnm="_receiver_black_";
73 > # Add receiver surfaces (rectangular)
74 > my $fmodnm="receiver_face";
75 > my $bmodnm="receiver_behind";
76   open(RADSCN, ">> $radscn");
77 < print RADSCN "void glow $modnm\n0\n0\n4 0 0 0 0\n\n";
78 < print RADSCN "$modnm polygon _receiver_\n0\n0\n12\n";
79 < print RADSCN "\t",$dim[0],"\t",$dim[2],"\t",$dim[5]+1e-5,"\n";
80 < print RADSCN "\t",$dim[0],"\t",$dim[3],"\t",$dim[5]+1e-5,"\n";
68 < print RADSCN "\t",$dim[1],"\t",$dim[3],"\t",$dim[5]+1e-5,"\n";
69 < print RADSCN "\t",$dim[1],"\t",$dim[2],"\t",$dim[5]+1e-5,"\n";
77 > print RADSCN "void glow $fmodnm\n0\n0\n4 1 1 1 0\n\n";
78 > print RADSCN "$fmodnm source f_receiver\n0\n0\n4 0 0 1 180\n";
79 > print RADSCN "void glow $bmodnm\n0\n0\n4 1 1 1 0\n\n";
80 > print RADSCN "$bmodnm source b_receiver\n0\n0\n4 0 0 -1 180\n";
81   close RADSCN;
82   # Generate octree
83   system "oconv -w $radscn > $octree";
84   die "Could not compile scene\n" if ( $? );
85 < # Set up sampling
85 > # Set up sampling of interior portal
86   # Kbin to produce incident direction in full Klems basis with (x1,x2) randoms
87   my $tcal = '
88   DEGREE : PI/180;
89 + sq(x) : x*x;
90   Kpola(r) : select(r+1, -5, 5, 15, 25, 35, 45, 55, 65, 75, 90);
91   Knaz(r) : select(r, 1, 8, 16, 20, 24, 24, 24, 16, 12);
92   Kaccum(r) : if(r-.5, Knaz(r) + Kaccum(r-1), 0);
# Line 85 | Line 97 | Kcol = Kbin - Kaccum(Krow-1);
97   Kazi = 360*DEGREE * (Kcol + (.5 - x2)) / Knaz(Krow);
98   Kpol = DEGREE * (x1*Kpola(Krow) + (1-x1)*Kpola(Krow-1));
99   sin_kpol = sin(Kpol);
100 < Dx = -cos(Kazi)*sin_kpol;
100 > Dx = cos(Kazi)*sin_kpol;
101   Dy = sin(Kazi)*sin_kpol;
102   Dz = sqrt(1 - sin_kpol*sin_kpol);
103 < Komega = 2*PI*if(Kbin-.5,
104 <        (cos(Kpola(Krow-1)*DEGREE) - cos(Kpola(Krow)*DEGREE))/Knaz(Krow),
105 <        1 - cos(Kpola(1)*DEGREE));
103 > KprojOmega = PI * if(Kbin-.5,
104 >        (sq(cos(Kpola(Krow-1)*DEGREE)) - sq(cos(Kpola(Krow)*DEGREE)))/Knaz(Krow),
105 >        1 - sq(cos(Kpola(1)*DEGREE)));
106   ';
107 < # Compute Klems bin from exiting ray direction
107 > # Compute Klems bin from exiting ray direction (forward or backward)
108   my $kcal = '
109   DEGREE : PI/180;
110 + abs(x) : if(x, x, -x);
111   Acos(x) : 1/DEGREE * if(x-1, 0, if(-1-x, 0, acos(x)));
112   posangle(a) : if(-a, a + 2*PI, a);
113   Atan2(y,x) : 1/DEGREE * posangle(atan2(y,x));
# Line 115 | Line 128 | kbin2(pol,azi) = select(kfindrow(1, pol),
128                  kaccum(7) + kazn(azi,360/knaz(8)),
129                  kaccum(8) + kazn(azi,360/knaz(9))
130          );
131 < kbin = kbin2(Acos(Dz), Atan2(Dy, -Dx));
131 > kbin = kbin2(Acos(abs(Dz)),Atan2(Dy,Dx));
132   ';
133   my $ndiv = 145;
134   my $nx = int(sqrt($nsamp*($dim[1]-$dim[0])/($dim[3]-$dim[2])) + .5);
135   my $ny = int($nsamp/$nx + .5);
136   $nsamp = $nx * $ny;
137   # Compute scattering data using rtcontrib
138 < my $cmd = "cnt $ndiv $ny $nx | rcalc -of -e '$tcal' " .
138 > my @tfarr;
139 > my @rfarr;
140 > my @tbarr;
141 > my @rbarr;
142 > my $cmd;
143 > my $rtcmd = "rtcontrib -h -ff -fo -n $nproc -c $nsamp " .
144 >        "-e '$kcal' -b kbin -bn $ndiv " .
145 >        "-o '$td/%s.flt' -m $fmodnm -m $bmodnm $rtargs $octree";
146 > my $rccmd = "rcalc -e '$tcal' " .
147 >        "-e 'mod(n,d):n-floor(n/d)*d' -e 'Kbin=mod(recno-.999,$ndiv)' " .
148 >        q{-if3 -e '$1=(0.265*$1+0.670*$2+0.065*$3)/KprojOmega'};
149 > if ( $doforw ) {
150 > $cmd = "cnt $ndiv $ny $nx | rcalc -of -e '$tcal' " .
151 >        "-e 'xp=(\$3+rand(.12*recno+288))*(($dim[1]-$dim[0])/$nx)+$dim[0]' " .
152 >        "-e 'yp=(\$2+rand(.37*recno-44))*(($dim[3]-$dim[2])/$ny)+$dim[2]' " .
153 >        "-e 'zp:$dim[4]' " .
154 >        q{-e 'Kbin=$1;x1=rand(2.75*recno+3.1);x2=rand(-2.01*recno-3.37)' } .
155 >        q{-e '$1=xp-Dx;$2=yp-Dy;$3=zp-Dz;$4=Dx;$5=Dy;$6=Dz' } .
156 >        "| $rtcmd";
157 > system "$cmd" || die "Failure running: $cmd\n";
158 > @tfarr = `$rccmd $td/$fmodnm.flt`;
159 > die "Failure running: $rccmd $td/$fmodnm.flt\n" if ( $? );
160 > @rfarr = `$rccmd $td/$bmodnm.flt`;
161 > die "Failure running: $rccmd $td/$bmodnm.flt\n" if ( $? );
162 > }
163 > if ( $doback ) {
164 > $cmd = "cnt $ndiv $ny $nx | rcalc -of -e '$tcal' " .
165          "-e 'xp=(\$3+rand(.35*recno-15))*(($dim[1]-$dim[0])/$nx)+$dim[0]' " .
166          "-e 'yp=(\$2+rand(.86*recno+11))*(($dim[3]-$dim[2])/$ny)+$dim[2]' " .
167 <        "-e 'zp:$dim[4]-1e-5' " .
167 >        "-e 'zp:$dim[5]' " .
168          q{-e 'Kbin=$1;x1=rand(1.21*recno+2.75);x2=rand(-3.55*recno-7.57)' } .
169 <        q{-e '$1=xp;$2=yp;$3=zp;$4=Dx;$5=Dy;$6=Dz' } .
170 <        "| rtcontrib -h -ff -n $nproc -c $nsamp -e '$kcal' -b kbin -bn $ndiv " .
171 <        "-m $modnm -w -ab 5 -ad 700 -lw 3e-6 $octree " .
172 <        "| rcalc -e 'x1:.5;x2:.5;$tcal' -e 'Kbin=floor((recno-1)/$ndiv)' " .
173 <        q{-if3 -e '$1=(0.265*$1+0.670*$2+0.065*$3)/(Komega*Dz)'};
174 < my @darr = `$cmd`;
175 < die "Failure running: $cmd\n" if ( $? );
169 >        q{-e '$1=xp-Dx;$2=yp-Dy;$3=zp+Dz;$4=Dx;$5=Dy;$6=-Dz' } .
170 >        "| $rtcmd";
171 > system "$cmd" || die "Failure running: $cmd\n";
172 > @tbarr = `$rccmd $td/$bmodnm.flt`;
173 > die "Failure running: $rccmd $td/$bmodnm.flt\n" if ( $? );
174 > @rbarr = `$rccmd $td/$fmodnm.flt`;
175 > die "Failure running: $rccmd $td/$fmodnm.flt\n" if ( $? );
176 > }
177   # Output XML prologue
178   print
179   '<?xml version="1.0" encoding="UTF-8"?>
# Line 236 | Line 276 | print '                        </Material>
276                          </AngleBasisBlock>
277                  </AngleBasis>
278          </DataDefinition>
279 <        <WavelengthData>
279 > ';
280 > if ( $doforw ) {
281 > print '         <WavelengthData>
282                  <LayerNumber>System</LayerNumber>
283                  <Wavelength unit="Integral">Visible</Wavelength>
284                  <SourceSpectrum>CIE Illuminant D65 1nm.ssp</SourceSpectrum>
# Line 248 | Line 290 | print '                        </Material>
290                          <ScatteringDataType>BTDF</ScatteringDataType>
291                          <ScatteringData>
292   ';
293 < # Output computed data (transposed order)
293 > # Output front transmission (transposed order)
294   for (my $od = 0; $od < $ndiv; $od++) {
295          for (my $id = 0; $id < $ndiv; $id++) {
296 <                print $darr[$ndiv*$id + $od];
296 >                print $tfarr[$ndiv*$id + $od];
297          }
298          print "\n";
299   }
258 # Output XML epilogue
300   print
301   '               </ScatteringData>
302          </WavelengthDataBlock>
303          </WavelengthData>
304 < </Layer>
304 >        <WavelengthData>
305 >                <LayerNumber>System</LayerNumber>
306 >                <Wavelength unit="Integral">Visible</Wavelength>
307 >                <SourceSpectrum>CIE Illuminant D65 1nm.ssp</SourceSpectrum>
308 >                <DetectorSpectrum>ASTM E308 1931 Y.dsp</DetectorSpectrum>
309 >                <WavelengthDataBlock>
310 >                        <WavelengthDataDirection>Reflection Front</WavelengthDataDirection>
311 >                        <ColumnAngleBasis>LBNL/Klems Full</ColumnAngleBasis>
312 >                        <RowAngleBasis>LBNL/Klems Full</RowAngleBasis>
313 >                        <ScatteringDataType>BRDF</ScatteringDataType>
314 >                        <ScatteringData>
315 > ';
316 > # Output front reflection (transposed order)
317 > for (my $od = 0; $od < $ndiv; $od++) {
318 >        for (my $id = 0; $id < $ndiv; $id++) {
319 >                print $rfarr[$ndiv*$id + $od];
320 >        }
321 >        print "\n";
322 > }
323 > print
324 > '               </ScatteringData>
325 >        </WavelengthDataBlock>
326 >        </WavelengthData>
327 > ';
328 > }
329 > if ( $doback ) {
330 > print '         <WavelengthData>
331 >                <LayerNumber>System</LayerNumber>
332 >                <Wavelength unit="Integral">Visible</Wavelength>
333 >                <SourceSpectrum>CIE Illuminant D65 1nm.ssp</SourceSpectrum>
334 >                <DetectorSpectrum>ASTM E308 1931 Y.dsp</DetectorSpectrum>
335 >                <WavelengthDataBlock>
336 >                        <WavelengthDataDirection>Transmission Back</WavelengthDataDirection>
337 >                        <ColumnAngleBasis>LBNL/Klems Full</ColumnAngleBasis>
338 >                        <RowAngleBasis>LBNL/Klems Full</RowAngleBasis>
339 >                        <ScatteringDataType>BTDF</ScatteringDataType>
340 >                        <ScatteringData>
341 > ';
342 > # Output back transmission (transposed order)
343 > for (my $od = 0; $od < $ndiv; $od++) {
344 >        for (my $id = 0; $id < $ndiv; $id++) {
345 >                print $tbarr[$ndiv*$id + $od];
346 >        }
347 >        print "\n";
348 > }
349 > print
350 > '               </ScatteringData>
351 >        </WavelengthDataBlock>
352 >        </WavelengthData>
353 >        <WavelengthData>
354 >                <LayerNumber>System</LayerNumber>
355 >                <Wavelength unit="Integral">Visible</Wavelength>
356 >                <SourceSpectrum>CIE Illuminant D65 1nm.ssp</SourceSpectrum>
357 >                <DetectorSpectrum>ASTM E308 1931 Y.dsp</DetectorSpectrum>
358 >                <WavelengthDataBlock>
359 >                        <WavelengthDataDirection>Reflection Back</WavelengthDataDirection>
360 >                        <ColumnAngleBasis>LBNL/Klems Full</ColumnAngleBasis>
361 >                        <RowAngleBasis>LBNL/Klems Full</RowAngleBasis>
362 >                        <ScatteringDataType>BRDF</ScatteringDataType>
363 >                        <ScatteringData>
364 > ';
365 > # Output back reflection (transposed order)
366 > for (my $od = 0; $od < $ndiv; $od++) {
367 >        for (my $id = 0; $id < $ndiv; $id++) {
368 >                print $rbarr[$ndiv*$id + $od];
369 >        }
370 >        print "\n";
371 > }
372 > print
373 > '               </ScatteringData>
374 >        </WavelengthDataBlock>
375 >        </WavelengthData>
376 > ';
377 > }
378 > # Output XML epilogue
379 > print '</Layer>
380   </Optical>
381   </WindowElement>
382   ';

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines