| 1 |
#!/bin/csh -f
|
| 2 |
# RCSid: $Id: vlpic.csh,v 3.4 2008/08/25 04:50:32 greg Exp $
|
| 3 |
#
|
| 4 |
# Compute falsecolor image of visibility level
|
| 5 |
# using the wacky formulas of Werner Adrian.
|
| 6 |
#
|
| 7 |
set age=40 # default age (years)
|
| 8 |
set tim=.2 # default time (seconds)
|
| 9 |
while ($#argv > 1)
|
| 10 |
switch ($argv[1])
|
| 11 |
case -a:
|
| 12 |
shift argv
|
| 13 |
set age="$argv[1]"
|
| 14 |
breaksw
|
| 15 |
case -t:
|
| 16 |
shift argv
|
| 17 |
set tim="$argv[1]"
|
| 18 |
breaksw
|
| 19 |
default:
|
| 20 |
echo bad option "'$argv[1]'"
|
| 21 |
exit 1
|
| 22 |
endsw
|
| 23 |
shift argv
|
| 24 |
end
|
| 25 |
set tc=`mktemp /tmp/vlcal.XXXXXXX`
|
| 26 |
set tp1=`mktemp /tmp/vlpic.XXXXXX`
|
| 27 |
set tp2=`mktemp /tmp/vlr2pic.XXXXXX`
|
| 28 |
set tp4=`mktemp /tmp/vlr4pic.XXXXXX`
|
| 29 |
set tp8=`mktemp /tmp/vlr8pic.XXXXXX`
|
| 30 |
set tf=($tc $tp1 $tp2 $tp4 $tp8)
|
| 31 |
set inpic=$argv[1]
|
| 32 |
onintr quit
|
| 33 |
set pr=`getinfo -d < $inpic`
|
| 34 |
set pr=($pr[4] $pr[2])
|
| 35 |
# ( vwright V < $inpic ; cat ) > $tc << _EOF_
|
| 36 |
cat > $tc << _EOF_
|
| 37 |
{ A : 3438 * sqrt(Vhn/xmax*Vvn/ymax); { pixel size (in minutes) } }
|
| 38 |
PI : 3.14159265358979323846;
|
| 39 |
A = 180*60/PI * sqrt(S(1)/PI);
|
| 40 |
sq(x) : x*x;
|
| 41 |
max(a,b) : if(a-b, a, b);
|
| 42 |
Lv : 0; { veiling luminance {temporarily zero} }
|
| 43 |
{ find direction of maximum contrast }
|
| 44 |
xoff(d) = select(d, 1, 1, 0, -1, -1, -1, 0, 1);
|
| 45 |
yoff(d) = select(d, 0, 1, 1, 1, 0, -1, -1, -1);
|
| 46 |
cf(lf, lb) = (lf - lb)/(lb + Lv);
|
| 47 |
contrast(d) = cf(li(1,xoff(d),yoff(d)), li(1,-xoff(d),-yoff(d)));
|
| 48 |
cwin(d1, d2) = if(contrast(d1) - contrast(d2), d1, d2);
|
| 49 |
bestdir(d) = if(d-1.5, cwin(d, bestdir(d-1)), d);
|
| 50 |
maxdir = bestdir(8);
|
| 51 |
Lt = WE*li(1,xoff(maxdir),yoff(maxdir)) + Lv;
|
| 52 |
La = WE*li(1,-xoff(maxdir),-yoff(maxdir)) + Lv;
|
| 53 |
LLa = log10(La);
|
| 54 |
F = if(La-.6, log10(4.1925) + LLa*.1556 + .1684*La^.5867,
|
| 55 |
if(La-.00418, 10^(-.072 + LLa*(.3372 + LLa*.0866)),
|
| 56 |
10^(.028 + .173*LLa)));
|
| 57 |
L = if(La-.6, .05946*La^.466,
|
| 58 |
if(La-.00418, 10^(-1.256 + .319*LLa),
|
| 59 |
10^(-.891 + LLa*(.5275 + LLa*.0227))));
|
| 60 |
LAt = log10(A) + .523;
|
| 61 |
AA = .36 - .0972*LAt*LAt/(LAt*(LAt - 2.513) + 2.7895);
|
| 62 |
LLat = LLa + 6;
|
| 63 |
AL = .355 - .1217*LLat*LLat/(LLat*(LLat - 10.4) + 52.28);
|
| 64 |
AZ = sqrt(AA*AA + AL*AL)/2.1;
|
| 65 |
DL1 = 2.6*sq(F/A + L);
|
| 66 |
M = if(La-.004, 10^-(10^-(if(La-.1,.125,.075)*sq(LLa+1) + .0245)), 0);
|
| 67 |
TGB = .6*La^-.1488;
|
| 68 |
FCP = 1 - M*A^-TGB/(2.4*(DL1*(AZ+2)/2));
|
| 69 |
DL2 = DL1*(AZ + T)/T;
|
| 70 |
FA = if(Age-64, sq(Age-56.6)/116.3 + 1.43, sq(Age-19)/2160 + .99);
|
| 71 |
DL3 = DL2 * FA;
|
| 72 |
DL4 = if(La-Lt, DL3*FCP, DL3);
|
| 73 |
lo = (Lt - La)/DL4; { Output VL }
|
| 74 |
_EOF_
|
| 75 |
pcomb -w -e "Age:$age;T:$tim" -f $tc -o $inpic > $tp1
|
| 76 |
pfilt -1 -x /2 -y /2 $inpic \
|
| 77 |
| pcomb -w -e "Age:$age;T:$tim" -f $tc -o - \
|
| 78 |
| pfilt -1 -r 1 -x $pr[1] -y $pr[2] > $tp2
|
| 79 |
pfilt -1 -x /4 -y /4 $inpic \
|
| 80 |
| pcomb -w -e "Age:$age;T:$tim" -f $tc -o - \
|
| 81 |
| pfilt -1 -r 1 -x $pr[1] -y $pr[2] > $tp4
|
| 82 |
pfilt -1 -x /8 -y /8 $inpic \
|
| 83 |
| pcomb -w -e "Age:$age;T:$tim" -f $tc -o - \
|
| 84 |
| pfilt -1 -r 1 -x $pr[1] -y $pr[2] > $tp8
|
| 85 |
pcomb -e 'max(a,b):if(a-b,a,b)' -e 'lo=max(gi(1),max(gi(2),max(gi(3),gi(4))))' \
|
| 86 |
$tp1 $tp2 $tp4 $tp8 \
|
| 87 |
| falsecolor -l VL -s 56 -m 1 -r 'if(v-1/8,1.6*8/7*(v-1/8)-.6,0)' \
|
| 88 |
-g 'if(v-1/8,if(v-.453125,1.6-1.6*8/7*(v-1/8),8/3*8/7*(v-1/8)),0)' \
|
| 89 |
-b 'if(v-1/8,1-8/3*8/7*(v-1/8),0)'
|
| 90 |
quit:
|
| 91 |
rm -f $tf
|
