:

# Create a vertical axis on pressure levels from a given axis
# description in meters.
# $1 input data file name
# $2 H_0 atmospheric scale height [default 8000m]

# On 29.04.19 10:39, Levke Caesar wrote:
# Der Druck in Aeolus bei einer Hoehe z ist
# p(z) = p_s * e^{-z/H_0}
# wobei H_0 = 8000m die atmospheric scale height ist.
# Der surface pressure p_s (bzw beim rausschreiben auch sfp) aendert
# sich aber von Zeitschritt zu Zeitschritt.
# -> here, wie simplify by set p_s fixed to 1000

# Example output:
#zaxistype = pressure
#size = 6
#units = mbar
#levels = 882 687 535 417 325 24
#lbounds = 1000 779 607 472 368 287
#ubounds = 779 607 472 368 287 2


# cdo zaxisdes $1 > zaxis-input.txt
# Sigh that write typically 2 Z-axis descriptions into 1 file,
# One for SurfaceVariable, and one for CellVariable.
# We must sort that out hand, for now :-(

makezaxis() {
    declare -a zlevels zlevelb
    zlevels=($1)
    zlevelb=($2)
    H_0="$3"
    p_s=1000

    plevels=''
    sep=""
    for z in ${zlevels[*]}
    do
	#plev=`echo $p_s '* e(-'$z'/'$H_0')' | bc -l`
	plev=`echo $z | awk '{printf("%1.0f", '$p_s' * exp(-$0/'$H_0'));}'`
	plevels="$plevels$sep$plev"
	sep=" "
    done
    #echo $plevels

    nbounds=${#zlevelb[@]}
    nboundsm1=$(($nbounds - 1))
    zlbounds="${zlevelb[@]:0:$nboundsm1}"
    #echo zlbounds $zlbounds
    zubounds="${zlevelb[*]:1:$nboundsm1}"
    #echo zubounds $zubounds

    plbounds=""
    sep=""
    for z in ${zlbounds[*]}
    do
	plb=`echo $z | awk '{printf("%1.0f", '$p_s' * exp(-$0/'$H_0'));}'`
	plbounds="$plbounds$sep$plb"
	sep=" "
    done
    pubounds=""
    sep=""
    for z in ${zubounds[*]}
    do
	pub=`echo $z | awk '{printf("%1.0f", '$p_s' * exp(-$0/'$H_0'));}'`
	pubounds="$pubounds$sep$pub"
	sep=" "
    done

echo 'zaxistype = pressure'
echo 'size = '${#zlevels[@]}
echo 'units = mbar'
echo 'levels = '${plevels}
echo 'lbounds = '${plbounds}
echo 'ubounds = '${pubounds}

}

declare -a levels_17 levelb_17
levels_17=(500 1500 2500 3500 4500 5500 6500 7500 8500 9500 10500 11500 12500 13500 14500 15500 33000)
levelb_17=(0 1000 2000 3000 4000 5000 6000 7000 8000 9000 10000 11000 12000 13000 14000 15000 16000 50000)

makezaxis "${levels_17[*]}" "${levelb_17[*]}" 8000 > zaxis_pressure_17

levels_6=(1000 3000 5000 7000 9000 30000)
levelb_6=(0 2000 4000 6000 8000 10000 50000)
makezaxis "${levels_6[*]}" "${levelb_6[*]}" 8000 > zaxis_pressure_6
