ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/util/genskyvec.pl
Revision: 2.2
Committed: Mon Jul 13 23:45:08 2009 UTC (13 years, 4 months ago) by greg
Content type: text/plain
Branch: MAIN
Changes since 2.1: +5 -4 lines
Log Message:
Fixed bug in ground sampling vectors

File Contents

# User Rev Content
1 greg 2.1 #!/usr/bin/perl -w
2 greg 2.2 # RCSid $Id: genskyvec.pl,v 2.1 2009/07/02 02:13:08 greg Exp $
3 greg 2.1 #
4     # Generate Reinhart vector for a given sky description
5     #
6     # G. Ward
7     #
8     use strict;
9     my @skycolor = (0.960, 1.004, 1.118);
10     while ($#ARGV >= 0) {
11     if ("$ARGV[0]" eq "-c") {
12     @skycolor = @ARGV[1..3];
13     shift @ARGV; shift @ARGV; shift @ARGV;
14     }
15     shift @ARGV;
16     }
17     # Load sky description into line array, separating sun if one
18     my @skydesc;
19     my $lightline;
20     my @sunval;
21     my $sunline;
22     my $skyOK = 0;
23     my $srcmod; # putting this inside loop breaks code(?!)
24     while (<>) {
25     push @skydesc, $_;
26     if (/^\w+\s+light\s+/) {
27     s/\s*$//; s/^.*\s//;
28     $srcmod = $_;
29     $lightline = $#skydesc;
30     } elsif (defined($srcmod) && /^($srcmod)\s+source\s/) {
31     @sunval = split(/\s+/, $skydesc[$lightline + 3]);
32     shift @sunval;
33     $sunline = $#skydesc;
34     } elsif (/\sskyfunc\s*$/) {
35     $skyOK = 1;
36     }
37     }
38     die "Bad sky description!\n" if (! $skyOK);
39     # Strip out the solar source if present
40     my @sundir;
41     if (defined $sunline) {
42     @sundir = split(/\s+/, $skydesc[$sunline + 3]);
43     shift @sundir;
44     undef @sundir if ($sundir[2] <= 0);
45     splice(@skydesc, $sunline, 5);
46     }
47     # Reinhart sky sample generator
48     my $rhcal = '
49     DEGREE : PI/180;
50     x1 = .5; x2 = .5;
51     MF : 2^2;
52     alpha : 90/(MF*7 + .5);
53     tnaz(r) : select(r, 30, 30, 24, 24, 18, 12, 6);
54     rnaz(r) : if(r-(7*MF-.5), 1, MF*tnaz(floor((r+.5)/MF) + 1));
55     raccum(r) : if(r-.5, rnaz(r-1) + raccum(r-1), 0);
56     RowMax : 7*MF + 1;
57     Rmax : raccum(RowMax);
58     Rfindrow(r, rem) : if(rem-rnaz(r)-.5, Rfindrow(r+1, rem-rnaz(r)), r);
59 greg 2.2 Rrow = if(Rbin-(Rmax-.5), RowMax-1, Rfindrow(0, Rbin));
60 greg 2.1 Rcol = Rbin - raccum(Rrow) - 1;
61     Razi_width = 2*PI / rnaz(Rrow);
62     RAH : alpha*DEGREE;
63     Razi = if(Rbin-.5, (Rcol + x2 - .5)*Razi_width, 2*PI*x2);
64 greg 2.2 Ralt = if(Rbin-.5, (Rrow + x1)*RAH, asin(-x1));
65 greg 2.1 Romega = if(.5-Rbin, 2*PI, if(Rmax-.5-Rbin,
66     Razi_width*(sin(RAH*(Rrow+1)) - sin(RAH*Rrow)),
67     2*PI*(1 - cos(RAH/2)) ) );
68     cos_ralt = cos(Ralt);
69     Dx = sin(Razi)*cos_ralt;
70     Dy = cos(Razi)*cos_ralt;
71     Dz = sin(Ralt);
72     ';
73     my $nbins = 2306; # This needs to be consistent with MF setting above
74     # Create octree for rtrace
75     my $octree = "/tmp/gtv$$.oct";
76     open OCONV, "| oconv - > $octree";
77     print OCONV @skydesc;
78     print OCONV "skyfunc glow skyglow 0 0 4 @skycolor 0\n";
79     print OCONV "skyglow source sky 0 0 4 0 0 1 360\n";
80     close OCONV;
81     # Run rtrace and average output for every 16 samples
82     my $tregcommand = "cnt $nbins 16 | rcalc -of -e '$rhcal' " .
83     q{-e 'Rbin=$1;x1=rand(recno*.37-5.3);x2=rand(recno*-1.47+.86)' } .
84     q{-e '$1=0;$2=0;$3=0;$4=Dx;$5=Dy;$6=Dz' } .
85     "| rtrace -h -ff -ab 0 -w $octree | total -if3 -16 -m";
86     my @tregval = `$tregcommand`;
87     unlink $octree;
88     # Find closest 3 patches to sun and divvy up direct solar contribution
89     my @bestdir;
90     if (@sundir) {
91     my $somega = ($sundir[3]/360)**2 * 3.141592654**3;
92 greg 2.2 my $cmd = "cnt " . ($nbins-1) .
93     " | rcalc -e '$rhcal' -e Rbin=recno " .
94 greg 2.1 "-e 'dot=Dx*$sundir[0] + Dy*$sundir[1] + Dz*$sundir[2]' " .
95     "-e 'cond=dot-.866' " .
96     q{-e '$1=if(1-dot,acos(dot),0);$2=Romega;$3=recno' };
97     @bestdir = `$cmd | sort -n | head -3`;
98     my (@ang, @dom, @ndx);
99     my $wtot = 0;
100     for my $i (0..2) {
101     ($ang[$i], $dom[$i], $ndx[$i]) = split(/\s+/, $bestdir[$i]);
102     $wtot += 1./($ang[$i]+.02);
103     }
104     for my $i (0..2) {
105     my $wt = 1./($ang[$i]+.02)/$wtot * $somega / $dom[$i];
106     my @scolor = split(/\s+/, $tregval[$ndx[$i]]);
107     for my $j (0..2) { $scolor[$j] += $wt * $sunval[$j]; }
108     $tregval[$ndx[$i]] = "@scolor\n";
109     }
110     }
111     # Output our final vector
112     print @tregval;