| 1 |
#!/bin/csh -f
|
| 2 |
# RCSid: $Id: phisteq.csh,v 3.3 2005/02/16 05:40:11 greg Exp $
|
| 3 |
set Ldmin=1 # minimum display luminance
|
| 4 |
set Ldmax=100 # maximum display luminance
|
| 5 |
set nsteps=100 # number of steps in perceptual histogram
|
| 6 |
set cvratio=0.08 # fraction of pixels to ignore in envelope clipping
|
| 7 |
set td=/tmp
|
| 8 |
set tf1=$td/hist$$
|
| 9 |
set tf1b=$td/hist$$.new
|
| 10 |
set tf2=$td/cumt$$
|
| 11 |
set tf3=$td/histeq$$.cal
|
| 12 |
set tf4=$td/cf$$.cal
|
| 13 |
set tf=($tf1 $tf1b $tf2 $tf3 $tf4)
|
| 14 |
if ( $#argv != 1 ) then
|
| 15 |
echo "Usage: $0 input.hdr > output.hdr"
|
| 16 |
exit 1
|
| 17 |
endif
|
| 18 |
set ifile=$1
|
| 19 |
set ibase=$ifile:t
|
| 20 |
if ( "$ibase" =~ *.pic ) set ibase=$ibase:r
|
| 21 |
set ibase=$ibase:t
|
| 22 |
onintr quit
|
| 23 |
cat > $tf3 << _EOF_
|
| 24 |
WE : 179; { Radiance white luminous efficacy }
|
| 25 |
Lmin : .0001; { minimum allowed luminance }
|
| 26 |
Ldmin : $Ldmin ; { minimum output luminance }
|
| 27 |
Ldmax : $Ldmax ; { maximum output luminance }
|
| 28 |
Stepsiz : 1/ $nsteps ; { brightness step size }
|
| 29 |
{ Daly local amplitude nonlinearity formulae }
|
| 30 |
sq(x) : x*x;
|
| 31 |
c1 : 12.6;
|
| 32 |
b : .63;
|
| 33 |
Bl(L) : L / (L + (c1*L)^b);
|
| 34 |
Lb(B) : (c1^b*B/(1-B))^(1/(1-b));
|
| 35 |
BLw(Lw) : Bl(Ldmin) + (Bl(Ldmax)-Bl(Ldmin))*cf(Bl(Lw));
|
| 36 |
{ first derivative functions }
|
| 37 |
Bl1(L) : (c1*L)^b*(1-b)/sq(L + (c1*L)^b);
|
| 38 |
Lb1(B) : c1^b/(1-b)/sq(1-B) * (c1^b*B/(1-B))^(b/(1-b));
|
| 39 |
{ derivative clamping function }
|
| 40 |
clamp2(L, aLw) : Lb(aLw) / L / Lb1(aLw) / (Bl(Ldmax)-Bl(Ldmin)) / Bl1(L);
|
| 41 |
clamp(L) : clamp2(L, BLw(L));
|
| 42 |
{ histogram equalization function }
|
| 43 |
lin = li(1);
|
| 44 |
Lw = WE/le(1) * lin;
|
| 45 |
Lout = Lb(BLw(Lw));
|
| 46 |
mult = if(Lw-Lmin, (Lout-Ldmin)/(Ldmax-Ldmin)/lin, 0) ;
|
| 47 |
ro = mult * ri(1);
|
| 48 |
go = mult * gi(1);
|
| 49 |
bo = mult * bi(1);
|
| 50 |
_EOF_
|
| 51 |
# Compute brightness histogram
|
| 52 |
pfilt -1 -p 1 -x 128 -y 128 $ifile | pvalue -o -b -d -h -H \
|
| 53 |
| rcalc -f $tf3 -e 'Lw=WE*$1;$1=if(Lw-Lmin,Bl(Lw),-1)' \
|
| 54 |
| histo 0 1 $nsteps | sed '/[ ]0$/d' > $tf1
|
| 55 |
# Clamp frequency distribution
|
| 56 |
set totcount=`sed 's/^.*[ ]//' $tf1 | total`
|
| 57 |
set tst=1
|
| 58 |
while ( $totcount > 0 )
|
| 59 |
sed 's/^.*[ ]//' $tf1 | total -1 -r \
|
| 60 |
| rcalc -e '$1=$1/'$totcount | rlam $tf1 - \
|
| 61 |
| tabfunc -i 0 cf > $tf4
|
| 62 |
if ( $tst <= 0 ) break
|
| 63 |
rcalc -f $tf4 -f $tf3 -e "T:$totcount*Stepsiz" \
|
| 64 |
-e 'clfq=floor(T*clamp(Lb($1)))' \
|
| 65 |
-e '$1=$1;$2=if($2-clfq,clfq,$2)' $tf1 > $tf1b
|
| 66 |
set newtot=`sed 's/^.*[ ]//' $tf1b | total`
|
| 67 |
set tst=`ev "floor((1-$cvratio)*$totcount)-$newtot"`
|
| 68 |
mv -f $tf1b $tf1
|
| 69 |
set totcount=$newtot
|
| 70 |
end
|
| 71 |
if ( $totcount < 1 ) then
|
| 72 |
# Fits in display range nicely already -- just normalize
|
| 73 |
pfilt -1 -e `pextrem $ifile | rcalc -e 'cond=recno-1.5;$1=1/(.265*$3+.67*$4+.065*$5)'` $ifile
|
| 74 |
else
|
| 75 |
# Plot the mapping function if we are in debug mode
|
| 76 |
if ( $?DEBUG ) then
|
| 77 |
cat > ${ibase}_histo.plt << _EOF_
|
| 78 |
include=curve.plt
|
| 79 |
title="Brightness Frequency Distribution"
|
| 80 |
subtitle= $ibase
|
| 81 |
ymin=0
|
| 82 |
xlabel="Perceptual Brightness B(Lw)"
|
| 83 |
ylabel="Frequency Count"
|
| 84 |
Alabel="Histogram"
|
| 85 |
Alintype=0
|
| 86 |
Blabel="Envelope"
|
| 87 |
Bsymsize=0
|
| 88 |
Adata=
|
| 89 |
_EOF_
|
| 90 |
(cat $tf1; echo \;; echo Bdata=) >> ${ibase}_histo.plt
|
| 91 |
rcalc -f $tf4 -f $tf3 -e "T:$totcount*Stepsiz" \
|
| 92 |
-e '$1=$1;$2=T*clamp(Lb($1))' $tf1 \
|
| 93 |
>> ${ibase}_histo.plt
|
| 94 |
cat > ${ibase}_brmap.plt << _EOF_
|
| 95 |
include=line.plt
|
| 96 |
title="Brightness Mapping Function"
|
| 97 |
subtitle= $ibase
|
| 98 |
xlabel="World Luminance (log cd/m^2)"
|
| 99 |
ylabel="Display Luminance (cd/m^2)"
|
| 100 |
ymax= $Ldmax
|
| 101 |
Adata=
|
| 102 |
_EOF_
|
| 103 |
cnt 100 | rcalc -f $tf4 -f $tf3 -e '$1=lx;$2=Lb(BLw(10^lx))' \
|
| 104 |
-e Lmin:Lb\(`sed -n '1s/[ ].*$//p' $tf1`\) \
|
| 105 |
-e Lmax:Lb\(`sed -n '$s/[ ].*$//p' $tf1`\) \
|
| 106 |
-e 'lx=$1/99*(log10(Lmax)-log10(Lmin))+log10(Lmin)' \
|
| 107 |
>> ${ibase}_brmap.plt
|
| 108 |
if ( $?DISPLAY ) then
|
| 109 |
bgraph ${ibase}_histo.plt ${ibase}_brmap.plt | x11meta &
|
| 110 |
endif
|
| 111 |
endif
|
| 112 |
# Map our picture
|
| 113 |
pcomb -f $tf4 -f $tf3 $ifile
|
| 114 |
endif
|
| 115 |
quit:
|
| 116 |
rm -f $tf
|
