#!/bin/sh
# script to plot 3-hour data from AWS for 1994
#
# each file contains 250 lines corresponding to measurements
# taken every 3 hours for the entire month
#
# algorithm: cat the files together then pipe through an
# awk script to add a decimal_day year  column and
# then send to xvgr to plot
#
# note: need to change day_start in awk script


echo "select a site to plot (q to quit):"
echo "kenton (k)"
read answer

case "$answer" in
	k) datafile=8922*
		parfile=kenton.par;;
	q) exit 1;;
esac

echo ""
echo "select data to plot: 10-minute (m) or 1-hour (h)? "
read answer  

case "$answer" in
	m) ;; # default 
	h) one_hour=1 ;; # set flag
esac

echo "plotting" "`basename $parfile .par`" "site ........."
echo ""

if [ -f "all.dat" ]
then
	echo "all.dat already exists; overwrite (y/n)? "
	read answer
fi

if [ ! -f "all.dat" -o "$answer" = y ]
then
	/usr/bin/rm all.dat

	# create all.dat 
		echo "creating file: all.dat"
		for i in ../raw/$datafile
		do
         gunzip -c --stdout $i | doit >> all.dat
         echo "gunzip -c --stdout" $i " | doit >> all.dat"
		done

fi

if [ -f "hr.dat" -a "$one_hour" -eq 1 ]
then
	echo "hr.dat already exists; overwrite (y/n)? "
	read answer
fi

if [ "$answer" = y ]
then
	/usr/bin/rm hr.dat

	# create hr.dat
	echo "newgreen all.dat > hr.dat"
	newgreen all.dat > hr.dat
	echo "cp hr.dat all.dat"
	cp hr.dat all.dat
fi

if [ -f "proc.dat" ]
then
	echo "output will be placed in file: proc.dat"
	echo "proc.dat already exists; overwrite (y/n)?"
	read answer
fi

if [ ! -f "proc.dat" -o "$answer" = y ]
then
	/usr/bin/rm proc.dat
	cat "all.dat" | \

	nawk 'BEGIN{ 
		a0=6.107799961;         b0=6.109177956
		a1=4.436518521e-1;      b1=5.034698970e-1
		a2=1.428945805e-2;      b2=1.886013408e-2
		a3=2.650648471e-4;      b3=4.176223716e-4
		a4=3.3031240396e-6;     b4=5.824720280e-6
		a5=2.034080948e-8;      b5=4.838803174e-8
		a6=6.136820929e-11;     b6=1.838826904e-10
		eps=0.622	
	} # end BEGIN

	{

		esat=a0+$2*(a1+$2*(a2+$2*(a3+$2*(a4+$2*(a5+a6*$2)))))   #sat Pv, mb
		e=esat*$6/100   # pressure of vapor, mb
		w=eps*e/($3-e)  # mixing ratio, mass vapor/mass dry air
		q=1000*w/(1+w)  # specific humidity, mass vapor/mass moist air (g/kg)
		
		printf("%4f %f %f %f %f %f %f %f\n", \
			$1, $2, $3, $4, $5, $6, q, $7)


	}' > proc.dat 
fi

# now plot it up

echo "Plot graph (y/n)?"
read answer

if [ "$answer" = y ]
then 

if [ -f "proc.dat" ]
then
	xvgr -p "$parfile" -source pipe -noredraw \
	-graph 0 \
		"nawk '{print "$"1, "$"2}' proc.dat" \
	-graph 1 \
		"nawk '{print "$"1, "$"3}' proc.dat" \
	-graph 2 \
		"nawk '{print "$"1, "$"4}' proc.dat" \
	-graph 3 \
		"nawk '{print "$"1, "$"5}' proc.dat" \
	-graph 4 \
		"nawk '{print "$"1, "$"6}' proc.dat" \
	-graph 5 \
		"nawk '{print "$"1, "$"7}' proc.dat" \
	-graph 6 \
		"nawk '{print "$"1, "$"8}' proc.dat" 

fi

fi