1 |
#!/usr/bin/perl -w |
2 |
# RCSid $Id: genklemsamp.pl,v 2.10 2011/04/16 00:39:07 greg Exp $ |
3 |
# |
4 |
# Sample Klems (full) directions impinging on surface(s) |
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 |
if ($#ARGV < 0) { |
12 |
print STDERR "Usage: genklemsamp [-c N ][-f{a|f|d}] [view opts] [geom.rad ..]\n"; |
13 |
exit 1; |
14 |
} |
15 |
my $nsamp = 1000; |
16 |
my $fmt = "a"; |
17 |
my @vopts="-vo 1e-5"; |
18 |
while ($#ARGV >= 0) { |
19 |
if ("$ARGV[0]" eq "-c") { |
20 |
shift @ARGV; |
21 |
$nsamp = $ARGV[0]; |
22 |
} elsif ("$ARGV[0]" =~ /^-f[afd]$/) { |
23 |
$fmt = substr($ARGV[0], 2, 1); |
24 |
} elsif ("$ARGV[0]" =~ /^-v[pdu]$/) { |
25 |
push @vopts, "@ARGV[0..3]"; |
26 |
shift @ARGV; shift @ARGV; shift @ARGV; |
27 |
} elsif ("$ARGV[0]" =~ /^-v[fhvo]$/) { |
28 |
push @vopts, "@ARGV[0..1]"; |
29 |
shift @ARGV; |
30 |
} elsif ("$ARGV[0]" =~ /^-v./) { |
31 |
print STDERR "Unsupported view option: $ARGV[0]\n"; |
32 |
exit 1; |
33 |
} elsif ("$ARGV[0]" =~ /^-./) { |
34 |
print STDERR "Unknown option: $ARGV[0]\n"; |
35 |
exit 1; |
36 |
} else { |
37 |
last; |
38 |
} |
39 |
shift @ARGV; |
40 |
} |
41 |
my $ofmt = ""; |
42 |
$ofmt = "-o$fmt" if ("$fmt" ne "a"); |
43 |
# Require parallel view and compile settings |
44 |
push @vopts, "-vtl"; |
45 |
my $vwset = `vwright @vopts V`; |
46 |
chomp $vwset; |
47 |
$_ = $vwset; s/^.*Vhn://; s/;.*$//; |
48 |
my $width = $_; |
49 |
$_ = $vwset; s/^.*Vvn://; s/;.*$//; |
50 |
my $height = $_; |
51 |
my $sca = sqrt($nsamp/($width*$height)); |
52 |
# Kbin is input to get direction in full Klems basis with (x1,x2) randoms |
53 |
my $tcal = 'DEGREE : PI/180;' . |
54 |
'Kpola(r) : select(r+1, -5, 5, 15, 25, 35, 45, 55, 65, 75, 90);' . |
55 |
'Knaz(r) : select(r, 1, 8, 16, 20, 24, 24, 24, 16, 12);' . |
56 |
'Kaccum(r) : if(r-.5, Knaz(r) + Kaccum(r-1), 0);' . |
57 |
'Kmax : Kaccum(Knaz(0));' . |
58 |
'Kfindrow(r, rem) : if(rem-Knaz(r)-.5, Kfindrow(r+1, rem-Knaz(r)), r);' . |
59 |
'Krow = if(Kbin-(Kmax+.5), 0, Kfindrow(1, Kbin));' . |
60 |
'Kcol = Kbin - Kaccum(Krow-1) - 1;' . |
61 |
'Kazi = 360*DEGREE * (Kcol + .5 - x2) / Knaz(Krow);' . |
62 |
'Kpol = DEGREE * (x1*Kpola(Krow) + (1-x1)*Kpola(Krow-1));' . |
63 |
'sin_kpol = sin(Kpol);' . |
64 |
'K0 = cos(Kazi)*sin_kpol;' . |
65 |
'K1 = sin(Kazi)*sin_kpol;' . |
66 |
'K2 = sqrt(1 - sin_kpol*sin_kpol);' ; |
67 |
my $ndiv = 145; |
68 |
my ($td, $inchk, $rmtmp, $orgcnt); |
69 |
if ($windoz) { |
70 |
my $tmploc = `echo \%TMP\%`; |
71 |
chomp($tmploc); |
72 |
$td = mkdtemp("$tmploc\\genklemsamp.XXXXXX"); |
73 |
chomp $td; |
74 |
$inchk = "rtrace -w -h -opLN $td\\surf.oct " . |
75 |
q{| rcalc -e "cond=and(5e9-$4,nOK);$1=$1;$2=$2;$3=$3" } . |
76 |
"-e \"and(a,b):if(a,b,a);sq(x):x*x\" -e \"$vwset\" " . |
77 |
q{-e "nOK=sq(Vdx*$5+Vdy*$6+Vdz*$7)-.999" } . |
78 |
" > $td\\origins.txt"; |
79 |
$rmtmp = "rmdir /S /Q $td"; |
80 |
} else { |
81 |
$td = mkdtemp("/tmp/genklemsamp.XXXXXX"); |
82 |
chomp $td; |
83 |
$inchk = "rtrace -w -h -faf -opLN $td/surf.oct " . |
84 |
q{| rcalc -e 'cond=and(5e9-$4,nOK);$1=$1;$2=$2;$3=$3' } . |
85 |
"-e 'and(a,b):if(a,b,a);sq(x):x*x' -e '$vwset' " . |
86 |
q{-e 'nOK=sq(Vdx*$5+Vdy*$6+Vdz*$7)-.999' } . |
87 |
"-if7 -of > $td/origins.flt"; |
88 |
$rmtmp = "rm -rf $td"; |
89 |
} |
90 |
sub GetNOrigins { |
91 |
if ($windoz) { |
92 |
my $count = 0; |
93 |
open(FILE, "< $td\\origins.txt"); |
94 |
$count++ while <FILE>; |
95 |
close(FILE); |
96 |
return $count; |
97 |
} |
98 |
return (-s "$td/origins.flt") / (4*3); |
99 |
} |
100 |
# Do we have any Radiance input files? |
101 |
if ($#ARGV >= 0) { |
102 |
system "oconv -f @ARGV > $td/surf.oct"; |
103 |
# Set our own view center and size based on bounding cube |
104 |
$_ = $vwset; s/^.*Vdx://; s/;.*$//; |
105 |
my @vd = $_; |
106 |
$_ = $vwset; s/^.*Vdy://; s/;.*$//; |
107 |
push @vd, $_; |
108 |
$_ = $vwset; s/^.*Vdz://; s/;.*$//; |
109 |
push @vd, $_; |
110 |
my @bcube = split ' ', `getinfo -d < $td/surf.oct`; |
111 |
$width = $bcube[3]*sqrt(3); |
112 |
$height = $width; |
113 |
push @vopts, ("-vp", $bcube[0]+$bcube[3]/2-$width/2*$vd[0], |
114 |
$bcube[1]+$bcube[3]/2-$width/2*$vd[1], |
115 |
$bcube[2]+$bcube[3]/2-$width/2*$vd[2]); |
116 |
push @vopts, ("-vh", $width, "-vv", $height); |
117 |
$vwset = `vwright @vopts V`; |
118 |
$sca = sqrt($nsamp/($width*$height)); |
119 |
my $xres; |
120 |
my $yres; |
121 |
my $ntot = 0; |
122 |
# This generally passes through the loop twice to get density right |
123 |
while ($ntot < $nsamp) { |
124 |
$xres = int($width*$sca) + 1; |
125 |
$yres = int($height*$sca) + 1; |
126 |
system "vwrays -x $xres -y $yres -pa 0 -pj .7 @vopts | $inchk"; |
127 |
$ntot = GetNOrigins(); |
128 |
if ($ntot == 0) { |
129 |
if ($nsamp < 200) { |
130 |
$sca = sqrt(200/($width*$height)); |
131 |
redo; |
132 |
} |
133 |
print STDERR "View direction does not correspond to any surfaces\n"; |
134 |
exit 1; |
135 |
} |
136 |
$sca *= 1.05 * sqrt($nsamp/$ntot); |
137 |
} |
138 |
# All set to produce our samples |
139 |
for (my $k = 1; $k <= $ndiv; $k++) { |
140 |
my $rn = rand(10); |
141 |
my $r1 = rand; my $r2 = rand; |
142 |
# Chance of using = (number_still_needed)/(number_left_avail) |
143 |
if ($windoz) { |
144 |
system "rcalc $ofmt -e \"$tcal\" -e \"$vwset\" " . |
145 |
"-e \"cond=($nsamp-outno+1)/($ntot-recno+1)-rand($rn+recno)\" " . |
146 |
"-e \"Kbin=$k;x1=rand($r1+recno);x2=rand($r2+recno)\" " . |
147 |
q{-e "$1=$1+Vo*Vdx; $2=$2+Vo*Vdy; $3=$3+Vo*Vdz" } . |
148 |
q{-e "$4=-K0*Vhx-K1*Vvx+K2*Vdx" } . |
149 |
q{-e "$5=-K0*Vhy-K1*Vvy+K2*Vdy" } . |
150 |
q{-e "$6=-K0*Vhz-K1*Vvz+K2*Vdz" } . |
151 |
"$td\\origins.txt"; |
152 |
} else { |
153 |
system "rcalc $ofmt -e '$tcal' -e '$vwset' " . |
154 |
"-e 'cond=($nsamp-outno+1)/($ntot-recno+1)-rand($rn+recno)' " . |
155 |
"-e 'Kbin=$k;x1=rand($r1+recno);x2=rand($r2+recno)' " . |
156 |
q{-e '$1=$1+Vo*Vdx; $2=$2+Vo*Vdy; $3=$3+Vo*Vdz' } . |
157 |
q{-e '$4=-K0*Vhx-K1*Vvx+K2*Vdx' } . |
158 |
q{-e '$5=-K0*Vhy-K1*Vvy+K2*Vdy' } . |
159 |
q{-e '$6=-K0*Vhz-K1*Vvz+K2*Vdz' } . |
160 |
"-if3 $td/origins.flt"; |
161 |
} |
162 |
} |
163 |
} else { |
164 |
# No Radiance input files, so just sample over parallel view |
165 |
my $xres = int($width*$sca) + 1; |
166 |
my $yres = int($height*$sca) + 1; |
167 |
my $ntot = $xres * $yres; |
168 |
for (my $k = 1; $k <= $ndiv; $k++) { |
169 |
my $rn = rand(10); |
170 |
my $r1 = rand; my $r2 = rand; |
171 |
my $r3 = rand; my $r4 = rand; |
172 |
# Chance of using = (number_still_needed)/(number_left_avail) |
173 |
if ($windoz) { |
174 |
system "cnt $yres $xres | rcalc $ofmt -e \"$tcal\" -e \"$vwset\" " . |
175 |
"-e \"xc=\$2;yc=\$1\" ". |
176 |
"-e \"cond=($nsamp-outno+1)/($ntot-recno+1)-rand($rn+recno)\" " . |
177 |
"-e \"Kbin=$k;x1=rand($r1+recno);x2=rand($r2+recno)\" " . |
178 |
"-e \"hpos=Vhn*((xc+rand($r3+recno))/$xres - .5)\" " . |
179 |
"-e \"vpos=Vvn*((yc+rand($r4+recno))/$yres - .5)\" " . |
180 |
q{-e "$1=Vpx+Vo*Vdx+hpos*Vhx+vpos*Vvx" } . |
181 |
q{-e "$2=Vpy+Vo*Vdy+hpos*Vhy+vpos*Vvy" } . |
182 |
q{-e "$3=Vpz+Vo*Vdz+hpos*Vhz+vpos*Vvz" } . |
183 |
q{-e "$4=-K0*Vhx-K1*Vvx+K2*Vdx" } . |
184 |
q{-e "$5=-K0*Vhy-K1*Vvy+K2*Vdy" } . |
185 |
q{-e "$6=-K0*Vhz-K1*Vvz+K2*Vdz" } ; |
186 |
} else { |
187 |
system "cnt $yres $xres | rcalc $ofmt -e '$tcal' -e '$vwset' " . |
188 |
q{-e 'xc=$2;yc=$1' } . |
189 |
"-e 'cond=($nsamp-outno+1)/($ntot-recno+1)-rand($rn+recno)' " . |
190 |
"-e 'Kbin=$k;x1=rand($r1+recno);x2=rand($r2+recno)' " . |
191 |
"-e 'hpos=Vhn*((xc+rand($r3+recno))/$xres - .5)' " . |
192 |
"-e 'vpos=Vvn*((yc+rand($r4+recno))/$yres - .5)' " . |
193 |
q{-e '$1=Vpx+Vo*Vdx+hpos*Vhx+vpos*Vvx' } . |
194 |
q{-e '$2=Vpy+Vo*Vdy+hpos*Vhy+vpos*Vvy' } . |
195 |
q{-e '$3=Vpz+Vo*Vdz+hpos*Vhz+vpos*Vvz' } . |
196 |
q{-e '$4=-K0*Vhx-K1*Vvx+K2*Vdx' } . |
197 |
q{-e '$5=-K0*Vhy-K1*Vvy+K2*Vdy' } . |
198 |
q{-e '$6=-K0*Vhz-K1*Vvz+K2*Vdz' } ; |
199 |
} |
200 |
} |
201 |
} |
202 |
system $rmtmp; |