#! /bin/bash
# Last edited on 2011-06-01 16:58:23 by stolfilocal

# Plots pressure core pressure against temperature, colored by time interva

# Arguments:
rix="$1"; shift # Reactor index (1, 2, or 3)
xtr="$1"; shift # Which temperature to use on X axis (0 = nozzle, 1 = bottom)
ypr="$1"; shift # Which pressure to use on Y axis (0 = coreA, 1 = coreB)

tma="$1"; shift # Tansition time A (hours since mar/11 00:00)
tmb="$1"; shift # Tansition time B (hours since mar/11 00:00)
tmc="$1"; shift # Tansition time C (hours since mar/11 00:00)
tmd="$1"; shift # Tansition time D (hours since mar/11 00:00)

# Internal variables:

utag="un${rix}"

tp_to_tag=( "TCn" "TCb" )
tp_to_tit=( "nozzle" "bottom" )

xtag="${tp_to_tag[${xtr}]}"
xtit="${tp_to_tit[${xtr}]}"
echo "${xtag}" "${xtit}"  1>&2

pr_to_tag=( "PCA" "PCB" )
pr_to_tit=( "core A" "core B" )

ytag="${pr_to_tag[${ypr}]}"
ytit="${pr_to_tit[${ypr}]}"
echo "${ytag}" "${ytit}" 1>&2

# Label style:

showlabs=0
labsty="left rotate by (-45) font \"arial,12\" offset 1,-0.7 tc rgb '#333333'"

if [[ showlabs -ne 0 ]]; then
  labplot="\"${tmp}.join\" using (dat(xc)):(dat(yc)):1 notitle with labels ${labsty},"
else
  labplot=""
fi

tmp="/tmp/$$"

pngfile="out/ptmp-${xtag}-${ytag}-${utag}-full.png"
  
echo "preparing the data file..." 1>&2

# Create a file with fields {time,PCA,PCB,TCn,TCb}: 
egrep -v -e '^[ ]*([\#]|$)' pres-${utag}-t.txt | sort -b -k5,5 > ${tmp}.data
egrep -v -e '^[ ]*([\#]|$)' heat-${utag}-t.txt | sort -b -k5,5 > ${tmp}.heat

join -1 5 -2 5 -o0,1.7,1.9,2.7,2.9 ${tmp}.data ${tmp}.heat \
  | gawk \
      -v xtr=${xtr} \
      -v ypr=${ypr} \
      '(($(2+ypr) == 99999) || ($(2+ypr) == 88888) || ($(4+xtr) == 99999) || ($(4+xtr) == 88888)) { next; } //{ print }' \
  | sort -b -k1,1g \
  > ${tmp}.join

export GDFONTPATH=..

echo "plotting..." 1>&2

gnuplot <<EOF
  set term png font arial 18 size 1660,1600
  set size ratio -1
  set output "${tmp}.png" 
  set xrange [-50:450]; unset logscale x; set xtics 50; set mxtics 5; set grid xtics
  set yrange [0.316:3162]; set logscale y; 
  set ytics ( \
    0.5, \
    1, \
    2.5, \
    5, \
    10, \
    25, \
    50, \
    100, \
    250, \
    500, \
    1000, \
    2500, \
    5000, \
    10000, \
    25000, \
    50000 \
  )
  set grid ytics
  set title "Fukushima Daiichi - Reactor #${rix}"
  set xlabel "${xtit} temperature ${xtag} (C)"
  set ylabel "absolute ${ytit} pressure ${ytag} (kPa abs)"
  set key left Left reverse
  set label "LIQUID" at  50,200 center font "arial,24"
  set label "VAPOR"  at 300,10  center font "arial,24"
  tc = 1;            # Column index of time coordinate.
  xc = 4 + ${xtr};   # Column index of X value (core pressure).
  yc = 2 + ${ypr};   # Column index of Y value (core temperature).
  tm0 = 0;
  tma = ${tma};
  tmb = ${tmb};
  tmc = ${tmc};
  tmd = ${tmd};
  tm9 = 100000;
  col(k) = column(k)
  dat(k) = ((col(k) == 99999) || (col(k) == 88888) ? 0/0 : col(k))
  din(k,t,a,b) = ((col(t) < a) || (col(t) > b) ? 0/0 : dat(k))
  plot \
    "steam-table.txt" using 1:3                       title "boiling"               with lines                  lc rgb '#ff7700', \
    (101.3)                                           title "atmos. press."         with lines                  lc rgb '#0077ff', \
    "${tmp}.join" using (dat(xc)):(dat(yc))           notitle                       with lines                  lc rgb '#444444', \
    ${labplot} \
    "${tmp}.join" using (dat(xc)):(din(yc,tc,tm0,tma)) title "0--${tma} hours"      with     points pt 7 ps 1.5 lc rgb '#cc0000', \
    "${tmp}.join" using (dat(xc)):(din(yc,tc,tma,tmb)) title "${tma}--${tmb} hours" with     points pt 7 ps 1.5 lc rgb '#aa8800', \
    "${tmp}.join" using (dat(xc)):(din(yc,tc,tmb,tmc)) title "${tmb}--${tmc} hours" with     points pt 7 ps 1.5 lc rgb '#339900', \
    "${tmp}.join" using (dat(xc)):(din(yc,tc,tmc,tmd)) title "${tmc}--${tmd} hours" with     points pt 7 ps 1.5 lc rgb '#0088ff', \
    "${tmp}.join" using (dat(xc)):(din(yc,tc,tmd,tm9)) title "after ${tmd} hours"   with     points pt 7 ps 1.5 lc rgb '#4400ff'
quit
EOF

echo "reducing..." 1>&2

convert ${tmp}.png -resize '50%' ${pngfile}

display ${pngfile}

rm -f ${tmp}.*

echo "done." 1>&2