#!/bin/csh -f

# sensitivity analysis using kuusk

# set shell variables to find codes
set ARCH = `uname -m`
set path = (. /home/plewis/bin/$ARCH  /home/plewis/bin/csh $path)

# some viewing/illumination parameters
set vza = 0
set sza = 45
set saa = 0

set name = (LAI soilBrightness clumping leafEccentricity NLeafLayers DryMatter senescenceProportion chlorophyll leafWater)

# check for required files
set required = (defaults.dat)
foreach i ($required)
  if ( ! -e $i ) then
    echo "required file $i" not found
    exit 0
  endif
end

# use: 
# -N number of leaf layers
# -l wavelength_start wavelength_end wavelength_step

# get data ranges
set lower = (0.001 0.05 0 -1  1 `gawk < defaults.dat '{print $4}'`)
set upper = (10    0.4  2  1  2.5 `gawk < defaults.dat '{print $5}'`)
set range = (`echo $lower $upper | gawk '{for(i=1;i<=NF/2;i++)print $(i+NF/2)-$i;}'`)

mkdir -p data
# define derivatives

foreach LAI (1.0)
foreach RSL1 (0.1)
foreach CLZ (0.0)
foreach EE (0.0)
foreach N (1.0)
foreach cdm (0.005)
foreach cs (0.1)
foreach cab (40)
foreach cw (0.012)

set params = ($LAI $RSL1 $CLZ $EE $N $cdm $cs $cab $cw)

set base = $params[1].$params[2].$params[3].$params[4].$params[5].$params[6].$params[7].$params[8].$params[9]
set test = `echo $params[4] |gawk '{if($1<0)print 0;else print 1;}'`
set testp = `echo $params[4] |gawk '{if($1<0)print -$1;else print $1;}'`
if ( $test == 0 ) then
    set ee = "-thm 90 -ee $testp"
else
    set ee = "-thm 0 -ee $testp"
endif
BRDF.ip -view $vza $vza 1 -solar_azimuth $sza -solar_zenith $sza -wavelength  400 2400 5 | kuusk -L $params[1] $ee -cAb $params[8] -rsl1 $params[2] -cW $params[9] -N $params[5] -clz $params[3] -cP $params[6] -cC $params[7] -data /home/plewis/src/kuusk/pigments.dat > tmp.$$
BRDF.ip -view $vza $vza 1  -solar_azimuth $sza -solar_zenith $sza -wavelength  400 2400 5 | kuusk -L $params[1] $ee -cAb $params[8] -rsl1 $params[2] -cW $params[9] -N $params[5] -clz $params[3] -cP $params[6] -cC $params[7] -data /home/plewis/src/kuusk/pigments.dat -ss > tmps.$$
BRDF.ip -view $vza $vza 1 -solar_azimuth $sza -solar_zenith $sza -wavelength  400 2400 5 | kuusk -L $params[1] $ee -cAb $params[8] -rsl1 $params[2] -cW $params[9] -N $params[5] -clz $params[3] -cP $params[6] -cC $params[7] -data /home/plewis/src/kuusk/pigments.dat -ms > tmpm.$$
head -1 < tmp.$$ | gawk '{for(i=3;i<=NF;i++)print $i}' > tmp2.$$
head -2 < tmp.$$ | tail -1 | gawk '{for(i=5;i<=NF;i++)print $i}' | paste tmp2.$$ - >  data/kuusk.$base.angles$vza.$sza.$saa
head -2 < tmps.$$ | tail -1 | gawk '{for(i=5;i<=NF;i++)print $i}' | paste tmp2.$$ - >  data/kuusk_sscat.$base.angles$vza.$sza.$saa
head -2 < tmpm.$$ | tail -1 | gawk '{for(i=5;i<=NF;i++)print $i}' | paste tmp2.$$ - >  data/kuusk_mscat.$base.angles$vza.$sza.$saa
\rm -f tmp*.$$ 

@ i = 1
while ( $i <= $#params )
  # derivative by fwd difference
  # reload
  set params = ($LAI $RSL1 $CLZ $EE $N $cdm $cs $cab $cw)
  set tmp = `echo $params[$i] $range[$i]| gawk '{print $1+$2*0.1}'`
  set dp = `echo $params[$i] $range[$i] | gawk '{print 0.1}'`
  if ( `echo $dp| gawk '{if($1==0)print 0;else print 1}'` == 1 ) then
    set params[$i] = $tmp
    # sort eccentricity
    set test = `echo $params[4] |gawk '{if($1<0)print 0;else print 1;}'`
    set testp = `echo $params[4] |gawk '{if($1<0)print -$1;else print $1;}'`
    if ( $test == 0 ) then
        set ee = "-thm 90 -ee $testp"
    else
	set ee = "-thm 0 -ee $testp"
    endif
    set fwd = $params[1].$params[2].$params[3].$params[4].$params[5].$params[6].$params[7].$params[8].$params[9]
    BRDF.ip -view $vza $vza 1 -solar_azimuth $sza -solar_zenith $sza -wavelength  400 2400 5 | kuusk -L $params[1] $ee -cAb $params[8] -rsl1 $params[2] -cW $params[9] -N $params[5] -clz $params[3] -cP $params[6] -cC $params[7] -data /home/plewis/src/kuusk/pigments.dat > tmp.$$
    head -1 < tmp.$$ | gawk '{for(i=3;i<=NF;i++)print $i}' > tmp2.$$
    BRDF.ip -view $vza $vza 1 -solar_azimuth $sza -solar_zenith $sza -wavelength  400 2400 5 | kuusk -L $params[1] $ee -cAb $params[8] -rsl1 $params[2] -cW $params[9] -N $params[5] -clz $params[3] -cP $params[6] -cC $params[7] -data /home/plewis/src/kuusk/pigments.dat -ss > tmps.$$
    BRDF.ip -view $vza $vza 1  -solar_azimuth $sza -solar_zenith $sza -wavelength  400 2400 5 | kuusk -L $params[1] $ee -cAb $params[8] -rsl1 $params[2] -cW $params[9] -N $params[5] -clz $params[3] -cP $params[6] -cC $params[7] -data /home/plewis/src/kuusk/pigments.dat -ms > tmpm.$$
    head -2 < tmp.$$  | tail -1| gawk '{for(i=5;i<=NF;i++)print $i}' | paste tmp2.$$ - >  data/kuusk.$fwd.angles$vza.$sza.$saa
    head -2 < tmps.$$  | tail -1| gawk '{for(i=5;i<=NF;i++)print $i}' | paste tmp2.$$ - >  data/kuusk_sscat.$fwd.angles$vza.$sza.$saa
    head -2 < tmpm.$$  | tail -1| gawk '{for(i=5;i<=NF;i++)print $i}' | paste tmp2.$$ - >  data/kuusk_mscat.$fwd.angles$vza.$sza.$saa
# do sscat/mscat
    \rm -f tmp2.$$ tmp.$$ tmps.$$ tmpm.$$
     paste  data/kuusk.$base.angles$vza.$sza.$saa  data/kuusk.$fwd.angles$vza.$sza.$saa | gawk '{print $1,($4-$2)/dp}' dp=$dp > data/kuuskderivative_parameter$i.$base.angles$vza.$sza.$saa
     paste  data/kuusk_sscat.$base.angles$vza.$sza.$saa  data/kuusk_sscat.$fwd.angles$vza.$sza.$saa | gawk '{print $1,($4-$2)/dp}' dp=$dp > data/kuuskderivative_sscat_parameter$i.$base.angles$vza.$sza.$saa

     paste  data/kuusk_mscat.$base.angles$vza.$sza.$saa  data/kuusk_mscat.$fwd.angles$vza.$sza.$saa | gawk '{print $1,($4-$2)/dp}' dp=$dp > data/kuuskderivative_mscat_parameter$i.$base.angles$vza.$sza.$saa
  endif

  @ i++
end
set params = ($LAI $RSL1 $CLZ $EE $N $cdm $cs $cab $cw)

# graph results
echo "base at $params"
@ i = 1
while ( $i <= $#params )
echo "parameter $i $name[$i]"
set params = ($LAI $RSL1 $CLZ $EE $N $cdm $cs $cab $cw)
set tmp = `echo $params[$i] $range[$i]| gawk '{print $1+$2*0.1}'`
set dp = `echo $params[$i] $range[$i] | gawk '{print $2*0.1}'`
if ( `echo $dp| gawk '{if($1==0)print 0;else print 1}'` == 1 ) then
    set params[$i] = $tmp
    set fwd = $params[1].$params[2].$params[3].$params[4].$params[5].$params[6].$params[7].$params[8].$params[9]
    echo "components"
    generate_graph   data/kuusk.$base.angles$vza.$sza.$saa data/kuusk.$fwd.angles$vza.$sza.$saa 
    generate_graph   data/kuusk_sscat.$base.angles$vza.$sza.$saa data/kuusk_sscat.$fwd.angles$vza.$sza.$saa data/kuusk_mscat.$base.angles$vza.$sza.$saa data/kuusk_mscat.$fwd.angles$vza.$sza.$saa
    echo "derivatives"
    generate_graph  data/kuuskderivative_sscat_parameter$i.$base.angles$vza.$sza.$saa data/kuuskderivative_mscat_parameter$i.$base.angles$vza.$sza.$saa  data/kuuskderivative_parameter$i.$base.angles$vza.$sza.$saa 
endif
@ i++
end

end
end
end
end
end
end
end
end
end