#!/bin/tcsh
#
# nplot [tab[3|6]_output_file [element_for distribution]]
#
# This script creates a GNUPLOT script for the plotting of EQ3NR scans,
# EQ6 titrations, or MINTEQA2 sweep runs (species distributions),
# and optionally also a PostScript version of the distribution graphs.
#
# These scans or titrations must have been processed by the "tab3" or "tab6"
# programs. They create a table with the actual value of the scanned variable
# (pH, Eh, T, p or a concentration) as 1st column, followed by the percentages
# of each species containing a given (user-defined) element.
# In case of EQ6 or MINTEQA2 output, also all precipitated solids
# are included.
#
# Depending on the first token of the first line in the distribution file,
# the x-axis scaling is equidistant or logarithmic.
#
# The name of the output files are the same as the corresponding input
# file, but instead of the extension ".[36M]d" they have ".[36M]gnu" for
# the gnuplot script file, and ".[36M]ps" for the PostScript graphics file.
# The output is optimized for gnuplot 3.7.
#
# The executable 'gnuplot' must be in the PATH.
#
# Author: Dr.Vinzenz Brendler, Forschungszentrum Dresden-Rossendorf e.V.
#                              Institute of Radiochemistry
#
# Version: 1.11                Date: 01-Sep-2010

onintr FINISH
stty icanon

if ($#argv < 1) then
	echo " "
	ls *.[36M]d >& /dev/null && ls *.[36M]d
	echo " "
	echo "Specify the name of the species distribution file: "
	echo " "
	echo -n "     < "
	set argv = (`gets`)
endif
if ($#argv > 2) then
	echo "Error: No more than two arguments are allowed: "
	echo "   species_distribution_file [ element ]"
	exit
endif
if (! -r $argv[1]) then
	echo "Error: Input file $argv[1] is not existent or not readable! "
	exit
endif
if (-z $argv[1]) then
	echo "Error: Input file $argv[1] is empty! "
	exit
endif

set CMD  = $0					# Command Version (nplot, nlogplot)
set DIST = $argv[1]				# Distribution File to be Processed
set TEMP = $DIST:r.$$			# Temporary File

# Determine whether the input was created by EQ3NR, EQ6 or MINTEQA2

set MARK = `sed -n '1s/^#\([36M]\).*/\1/p' $DIST`
switch ($MARK)
case "3":
	set CODE = EQ3NR; breaksw
case "6": 
	set CODE = EQ6; breaksw
case "M": 
	set CODE = MINTEQA2; breaksw
default:
	echo "Originating speciation code unknown ! " ; exit; breaksw
endsw

set SCRIPT = $DIST:r.${MARK}gnu	# Gnuplot Script File
set POST   = $DIST:r.${MARK}ps	# PostScript Graphic File

# Check whether the output file (Gnuplot Script) already exists

if (-e $SCRIPT) then
	echo -n "Gnuplot Script file $SCRIPT already exists, overwrite y/(n) ? "
	set rr = (`gets`)
	if ($rr != y && $rr != Y) exit
endif

# Ask for the element wanted for distribution

set SYMBOL = U				# Uranum as Default Species for Distribution

if ($#argv > 1) then
	set ELEMENT = $argv[2]
else
	echo " "
	echo -n "Element wanted for distribution ( <CR> = $SYMBOL ) : "
	set ELEMENT = (`gets`)
	if ( ! $#ELEMENT) set ELEMENT = "$SYMBOL"
endif

# Ask for the minimum fraction of concentration to be plotted

set DEFFRAC = 1.0	 			# Default Minimum Fraction (%) to be plotted

echo -n "Minimum speciation fraction (in %) for plotting ( <CR> = $DEFFRAC ) : "
set MINFRAC = (`gets`)
if ( ! $#MINFRAC) set MINFRAC = $DEFFRAC

# Ask for the plot title

echo -n "Title of the plot : "
set TITLE = (`gets`)
if ( ! $#TITLE) set TITLE = "Species distribution from $CODE"
echo " "

# Build up a list with Column# and SpeciesName for each species giving
# a speciation fraction larger than the user-defined minimum

awk 'BEGIN   { LIMIT = '$MINFRAC' } \
	NR == 1 { COLS = NF ; \
			for (i = 2; i <= COLS; i++) \
				{ name[i] = $i ; maxi[i] = 0.0 } ; next } \
		   { for (i = 2; i <= COLS; i++) \
				if ($i + 0 > maxi[i]) maxi[i] = $i } \
	END	   { for (i = 2; i <= COLS; i++) \
				if (maxi[i] + 0 > LIMIT) printf "%2d %s\n", i, name[i] } \
' $DIST > $TEMP

# Determine number of species to be plotted, and the scanned variable

set NSPECS = `cat $TEMP | wc -l`
set SCAN   = `sed -n '1s/#[36MEL]*#\([^ 	]*\)[ 	].*/\1/p' $DIST`

if ($NSPECS < 1) then
	echo "There are no species to plot"
	exit
endif

# Check for logarithmic plot request

set LOG = `sed -n '1s/^#[36M]\(.\).*/\1/p' $DIST`

# Build up the GNUPLOT script

echo "set title" \'$TITLE\'                 > $SCRIPT
echo "# Computed with" $CODE                >> $SCRIPT
echo "# Concentration cutoff at" $MINFRAC % >> $SCRIPT
echo "set xlabel" \'$SCAN\'                 >> $SCRIPT
echo "set ylabel" \'% $ELEMENT\'            >> $SCRIPT
echo "set yrange [0:100]"                   >> $SCRIPT
echo "set ytics 0,10,100"                   >> $SCRIPT
echo "set data style lines"                 >> $SCRIPT
if ($LOG == "L") then
     echo "set logscale x"                  >> $SCRIPT
endif

set i = 1
echo "plot" \\ >> $SCRIPT
while ($i <= $NSPECS)
	set COL     = `awk 'NR == '$i' {print $1}' $TEMP`
	set SPECIES = `awk 'NR == '$i' {print $2}' $TEMP`
	echo -n  \'$DIST\' us 1:$COL title \'$SPECIES\' >> $SCRIPT
	echo -n "."
	if ($i < $NSPECS) echo , \\ >> $SCRIPT
	@ i++
end
echo " " >> $SCRIPT
echo "pause -1" \'Hit return to continue\' >> $SCRIPT
echo " "

# Ask whether also a PostScript file should be created

# echo " "
# echo -n "Do you also want to create a PostScript file y/(n) ? "
# set rr = (`gets`)
# if ($rr == y || $rr == Y) set PS_OUT

# Append PostScript part to gnuplot script

if ($?PS_OUT) then
	echo "set terminal postscript landscape color" >> $SCRIPT
	echo "set output" \'$POST\' >> $SCRIPT
	echo "replot" >> $SCRIPT
endif
echo "        Successfully created: $SCRIPT"

# Run the gnuplot script and show graphs
# (skipped for alphanumeric terminals)

echo " "
if ($?DISPLAY) then
     if ($?PS_OUT) then
          gnuplot $SCRIPT && echo "        Successfully created: $POST"
          echo " "
          echo -n "Do you want to view the PostScript file y/(n) ? "
          set rr = (`gets`)
          if ($rr == y || $rr == Y) ghostview -a4 -landscape -magstep -2 $POST
     else
          gnuplot $SCRIPT
     endif
else
     echo "GNUPLOT can not be run, DISPLAY variable is not set."
endif

# Delete all temporary files

FINISH:
/bin/rm -f $TEMP

# Future Enhancements:
#