#!/bin/bash
### dclrinexplot.sh
### run fortran program dclrinex and produces plots with GMT
### version 20191219


### functions

_usage()
{
    cat << END
Usage:   dclrinex.sh <sss1> <sss2> <yydoy[.xxxx]> <dur> "[<option>]"

 where:

           sss1 = 4- or 9-character name for system 1 (travel)
           sss2 = 4- or 9-character name for system 2 (local)
          yydoy = first day of Rinex
          .xxxx = option to specify precise start time
            dur = Number of days to process
         option = lzXXXX to use LZ files with name lzXXXXmj.day for system 2
                  P1ISC1 to use C1 data as P1 when P1 is not available
                  default is none of the two options
         More options available when running the dclrinex program by itself
    Output file = sss1sss2yydoy_dur.ext
 Warning: Remove files which have the same filename as this output file!
END
}

_tmp_cleanup()
{
    rm -f *[CcPpEeBbRrDd][1235c].tmp
}

_get_lin_axes() # $FNAME
{
    FNAME=${1}

    ymax=$(awk -v tmax=-99999999 '{ if ($2 > tmax) tmax=$2 } END { print tmax } ' $FNAME)
    ymin=$(awk -v tmin=99999999  '{ if ($2 < tmin) tmin=$2 } END { print tmin } ' $FNAME)
    ymax=$(echo $ymax | awk '{ print int($1)+1 }')
    ymin=$(echo $ymin | awk '{ print int($1)-1 }')

    yspan=$(echo $ymax $ymin | awk '{ print $1-$2 }')

    yint=$(echo $yspan | awk '{ print $1/30 }')

    yatic=100
    if [[ ${yspan} -lt 200 ]]; then
        yatic=10
    fi
    if [[ ${yspan} -lt 20 ]]; then
        yatic=1
    fi
    yftic=$(echo ${yatic} | awk '{ print $1/5 }')

    xmax=$(awk -v tmax=-99999999 ' { if ($1 > tmax) tmax=$1 } END { print tmax } ' $FNAME)
    xmin=$(awk -v tmin=99999999  ' { if ($1 < tmin) tmin=$1 } END { print tmin } ' $FNAME)
    xmax=$(echo ${xmax} | awk '{ print int($1)+1 }')
    xmin=$(echo ${xmin} | awk '{ print int($1) }')

    xspan=$(echo ${xmax} ${xmin} | awk '{ print $1-$2 }')
    xftic=$(echo ${xspan} | awk '{ print int($1/10)+1 }')
    xatic=$(echo ${xftic} | awk '{ print $1*2 }')
    xint=$(echo ${xspan} | awk '{ print $1/30 }')
}

_get_log_axes() # $FNAME
{
    FNAME=${1}

    xmax=$(awk -v tmax=0 '{ if ($1 > tmax) tmax=$1 } END { print tmax } ' ${FNAME})
    xmin=$(awk -v tmin=99999999 ' { if ($1 < tmin) tmin=$1 } END { print tmin } ' ${FNAME})
    ymax=$(awk -v tmax=0 ' { if ($2 > tmax) tmax=$2 } END { print tmax } ' ${FNAME})
    ymin=$(awk -v tmin=99999999 ' { if ($2 < tmin) tmin=$2 } END { print tmin } ' ${FNAME})

    xmax=$(echo ${xmax}  | awk '{ print 10^(int(log($1)/2.302585093)+1) }' )
    xmin=$(echo ${xmin}  | awk '{ print 10^(int(log($1)/2.302585093)) }' )
    ymax=$(echo ${ymax}  | awk '{ print 10^(int(log($1)/2.302585093)) }' )
    ymin=$(echo ${ymin}  | awk '{ print 10^(int(log($1)/2.302585093)-1) }' )

    xspan=$(echo ${xmax} ${xmin} | awk '{ print $1/$2 }' )
    xftic=$(echo ${xmax} ${xspan} | awk '{ print $1/$2 }' )
    xatic=$(echo ${xftic}      | awk '{ print $1/1  }' )
    xgtic=$(echo ${xftic}      | awk '{ print $1/10  }' )

    yspan=$(echo ${ymax} ${ymin} | awk '{ print $1/$2 }' )
    yftic=$(echo ${ymax} ${yspan} | awk '{ print $1/$2 }' )
    yatic=$(echo ${yftic}      | awk '{ print $1/1  }' )
    ygtic=$(echo ${yftic}      | awk '{ print $1/10  }' )
}

_get_med_text() # $FNAME $SIZE $POS $DEVS_LIS
{
    FNAME=${1}
    SIZE=${2}
    POS=${3}
    DEVS_LIS=${4}

    cat ${FNAME} |
    while read x
    do
        med=$(echo ${x} | awk '{ print $3 }')
        Med=$(echo ${med} | awk '{ printf "%7.2f", $1 }')
        xpos=$(echo ${xmin} ${xint} | awk '{ print $1+2*$2 }')
        if [[ "x${POS}" = "xBOT" ]]; then
            ypos=$(echo ${ymin} ${ymax} ${yint} | awk '{ print 1.*$1+0*$2+1*$3}')
        elif [[ "x${POS}" = "xTOP" ]]; then
            ypos=$(echo ${ymin} ${ymax} ${yint} | awk '{ print 0.*$1+1*$2-1*$3}')
        fi
        if [[ -z ${cmd} ]]; then
            echo ${xpos} ${ypos} ${SIZE} 0 1 BL Median = ${Med} ns >> ${DEVS_LIS}
        else
            echo ${xpos} ${ypos} ${SIZE},1,blue BL Median = ${Med} ns >> ${DEVS_LIS}
        fi
    done
}

_get_tdev_text() # $FNAME $SIZE $CODE $POS $DEVS_LIS
{
    FNAME=${1}
    SIZE=${2}
    CODE=${3}
    POS=${4}
    DEVS_LIS=${5}

    i=0
    cat ${FNAME} |
    while read x
    do
        tau=$(echo ${x} | awk '{ print $1 }')
        tdev=$(echo ${x} | awk '{ print $2 }')
        Tau=$(echo ${tau} | awk '{ printf "%6.0f", $1 }')
        Tdev=$(echo ${tdev} | awk '{ printf "%6.0f", $1*1.e12 }')
        if [[ "x$POS" = "xLEFT" ]]; then
            xpos=$(echo ${xmax} ${xint} | awk '{ print $1+$2*5 }' )
        elif [[ "x$POS" = "xRIGHT" ]]; then
            xpos=$(echo ${xmax} ${xint} | awk '{ print $1+$2*12 }' )
        fi
        ypos=$(echo ${ymax} ${yint} ${i} | awk '{ printf "%9.2f", $1+$2*($3-4) }' )
        if [[ -z ${cmd} ]]; then
            echo ${xpos} ${ypos} ${SIZE} 0 1 LM ${Tau} s: ${CODE}= ${Tdev} ps >> ${DEVS_LIS}
        else
            echo ${xpos} ${ypos} ${SIZE},1,blue LM ${Tau} s: ${CODE}= ${Tdev} ps >> ${DEVS_LIS}
        fi
        i=$((${i} + 1))
    done
}

_gmt4_clear_all()
{
    rm .gmt* 2> /dev/null
    rm gmt.conf 2> /dev/null
    rm gmt.history 2> /dev/null
}

_gmt4_set()
{
    gmtset N_HEADER_RECS 0
    gmtset BASEMAP_AXES WeSn
    gmtset LABEL_FONT_SIZE 16p
    gmtset HEADER_FONT_SIZE 20p
    gmtset PAPER_MEDIA a4
}

_gmt5_clear_all()
{
    ${cmd} clear history 2>/dev/null
    ${cmd} clear conf 2>/dev/null
    ${cmd} clear cache 2>/dev/null
}

_gmt5_set()
{
    ${cmd} gmtset IO_HEADER 0
    ${cmd} gmtset MAP_FRAME_AXES WeSn
    ${cmd} gmtset FONT_LABEL 16p
    ${cmd} gmtset FONT_TITLE 20p
    ${cmd} gmtset PS_MEDIA a4
}

_plot_code() # $COD1 $COD2
{
    COD1=${1}
    COD2=${2}

    ### diff plot
    if [[ ! -z ${COD2} ]]; then bas=""; else bas="MJD"; fi

    # COD1 diff plot
    _get_lin_axes ${COD1}.tmp
    ${cmd} psxy ${COD1}.tmp -S-0.04 -Xa1.00i -Ya4.00i -R${xmin}/${xmax}/${ymin}/${ymax} -JX6.10i/2.50i -Bf${xftic}a${xatic}:"${bas}":/f${yftic}a${yatic}:"${COD1} /ns"::."${TITLE}": -W1,0/0/0 -K > ${OUTFILE}.ps

    # COD2 diff plot
    if [[ ! -z ${COD2} ]]; then
        _get_lin_axes ${COD2}.tmp
        ${cmd} psxy -O ${COD2}.tmp -S-0.04 -Xa1.00i -Ya1.00i -R${xmin}/${xmax}/${ymin}/${ymax} -JX6.10i/2.50i -Bf${xftic}a${xatic}:"MJD":/f${yftic}a${yatic}:"${COD2} /ns": -W1,0/0/0 -K >> ${OUTFILE}.ps
    fi


    ### Add text
    rm -f devs.lis
    # med
    _get_med_text med${COD1,,}.tmp 12 TOP devs.lis
    if [[ ! -z ${COD2} ]]; then
        _get_med_text med${COD2,,}.tmp 12 BOT devs.lis
    fi
    # tdev
    _get_tdev_text tdev${COD1,,}.tmp 9 ${COD1} LEFT  devs.lis
    if [[ ! -z ${COD2} ]]; then
        _get_tdev_text tdev${COD2,,}.tmp 9 ${COD2} RIGHT devs.lis
    fi
    #
    if [[ ! -z ${cmd} ]]; then opt="-F+f+j"; else opt="-G0/0/255"; fi
    ${cmd} pstext -N -O -K -R${xmin}/${xmax}/${ymin}/${ymax} -Xa1.00i -Ya1.00i -JX6.20i/5.50i ${opt} devs.lis >> ${OUTFILE}.ps
    rm -f devs.lis

    ### TDEV plot
    if [[ ! -z ${COD2} ]]; then bas=""; kon="-K"; else bas="Tau /s"; kon=""; fi

    # COD1 TDEV plot
    _get_log_axes tdev${COD1,,}.tmp
    ${cmd} psxy tdev${COD1,,}.tmp -O ${kon} -Xa8.50i -Ya3.70i -Sc0.05i -R${xmin}/${xmax}/${ymin}/${ymax} -JX2.50il/1.75il -Bf${xftic}a${xatic}:"${bas}":/f1a1:"Tdev${COD1} /s": -W5,0/0/0 >> ${OUTFILE}.ps

    # COD2 TDEV plot
    if [[ ! -z ${COD2} ]]; then
        _get_log_axes tdev${COD2,,}.tmp
        ${cmd} psxy tdev${COD2,,}.tmp -O -Xa8.50i -Ya1.00i -Sc0.05i -R${xmin}/${xmax}/${ymin}/${ymax} -JX2.50il/1.75il -Bf${xftic}a${xatic}:"Tau /s":/f1a1:"Tdev${COD2} /s": -W5,0/0/0 >> ${OUTFILE}.ps
    fi

    # transform to pdfr
    ps2pdf ${OUTFILE}.ps
    rm -f ${OUTFILE}.ps
    mv ${OUTFILE}.pdf ${OUTFILE}_${COD1}${COD2:1:1}.pdf

}


### arguments
if [[ $# -lt 4 ]]; then
    _usage
    exit 1
fi

### variables
doyd=$(echo ${3} | cut -c1-5)
OUTFILE=$(echo ${1:0:4}${2:0:4}${doyd}"_"${4})
DATERUN=$(date '+%F')
TITLE=$(echo ${DATERUN} ${OUTFILE})

# check GMT software version 4/5
cmd=$(command -v gmt)

if [[ -z ${cmd} ]]; then
    _gmt4_clear_all
    _gmt4_set
else
    _gmt5_clear_all
    _gmt5_set
fi


### core

# Delete old *.tmp files
_tmp_cleanup


# Invoke dclrinex
dclrinex "${1}" "${2}" "${3}" "${4}" "${5}"

# GPS
# plot P1/P2
if [[ -s "P1.tmp" && -s "P2.tmp" ]]; then
    _plot_code P1 P2
elif [[ -s "P1.tmp" && ! -s "P2.tmp" ]]; then
    _plot_code P1
elif [[ ! -s "P1.tmp" && -s "P2.tmp" ]]; then
    _plot_code P2
fi

# plot C1/C2
if [[ -s "C1.tmp" && -s "C2.tmp" ]]; then
    _plot_code C1 C2
elif [[ -s "C1.tmp" && ! -s "C2.tmp" ]]; then
    _plot_code C1
elif [[ ! -s "C1.tmp" && -s "C2.tmp" ]]; then
    _plot_code C2
fi

# plot C5
if [[ -s "C5.tmp" ]]; then
    _plot_code C5
fi

# GAL
# plot E1/E5
if [[ -s "E1.tmp" && -s "E5.tmp" ]]; then
    _plot_code E1 E5
elif [[ -s "E1.tmp" && ! -s "E5.tmp" ]]; then
    _plot_code E1
elif [[ ! -s "E1.tmp" && -s "E5.tmp" ]]; then
    _plot_code E5
fi

# BDS
# plot B1/B2
if [[ -s "B1.tmp" && -s "B2.tmp" ]]; then
    _plot_code B1 B2
elif [[ -s "B1.tmp" && ! -s "B2.tmp" ]]; then
    _plot_code B1
elif [[ ! -s "B1.tmp" && -s "B2.tmp" ]]; then
    _plot_code B2
fi

# plot B3/B5
if [[ -s "B3.tmp" && -s "B5.tmp" ]]; then
    _plot_code B3 B5
elif [[ -s "B3.tmp" && ! -s "B5.tmp" ]]; then
    _plot_code B3
elif [[ ! -s "B3.tmp" && -s "B5.tmp" ]]; then
    _plot_code B5
fi

# GLO P1/P2
# plot R1/R2
if [[ -s "R1.tmp" && -s "R2.tmp" ]]; then
    _plot_code R1 R2
elif [[ -s "R1.tmp" && ! -s "R2.tmp" ]]; then
    _plot_code R1
elif [[ ! -s "R1.tmp" && -s "R2.tmp" ]]; then
    _plot_code R2
fi

# GLO C1/C2
# plot D1/D2
if [[ -s "D1.tmp" && -s "D2.tmp" ]]; then
    _plot_code D1 D2
elif [[ -s "D1.tmp" && ! -s "D2.tmp" ]]; then
    _plot_code D1
elif [[ ! -s "D1.tmp" && -s "D2.tmp" ]]; then
    _plot_code D2
fi

echo "Output files: "${OUTFILE}

# clear GMT files
if [[ -z ${cmd} ]]; then
    _gmt4_clear_all
else
    _gmt5_clear_all
fi

### EOS
exit 0
