#!/bin/bash

	PARFILE=$1
	TRAXIS=$2
    MASKOUT=$3
    source $PARFILE
    
    ln -s /radar3/flex/age/age_mask/${TRAXIS} traxis.xy
    
 	BOUNDS="${WLON}/${ELON}/${SLAT}/${NLAT}"
 	
 	#	SWITCH GMT VERSION MANUALLY
	gmtswitch /Users/esg006/gmt5    
    
    gmtset PROJ_ELLIPSOID = Sphere
    
	ln -s ../../age/age.3.6.nc age.nc
	
	gmt grdmath age.nc 100 DIV = age_scaled.grd
	
	gmt grdsample age_scaled.grd -R$BOUNDS -I1m/1m -r -nl -Gage_tmp1.grd
	
	gmt grdinfo age.grd -C | awk '{print $2-0.1,$3+0.1,$4-0.1,$5+0.1}' > wesn.txt
	read WLONA ELONA SLATA NLATA < wesn.txt 
	
	ABOUNDS="${WLONA}/${ELONA}/${SLATA}/${NLATA}"
	
	gmt gmtselect traxis.xy -R$ABOUNDS > age_mask.xy 
	
	case "$MASKOUT" in
	
	N)	echo "Masking out to the NORTH of the trench"
		echo "${ELONA} ${NLATA}" >> age_mask.xy
		echo "${WLONA} ${NLATA}" >> age_mask.xy
		head -n1 age_mask.xy >> age_mask.xy
		;;
	E)	echo "Masking out to the EAST of the trench"
		echo "${ELONA} ${SLATA}" >> age_mask.xy
		echo "${ELONA} ${NLATA}" >> age_mask.xy
		head -n1 age_mask.xy >> age_mask.xy
		;;
	S)	echo "Masking out to the SOUTH of the trench"
		echo "${ELONA} ${SLATA}" >> age_mask.xy
		echo "${WLONA} ${SLATA}" >> age_mask.xy
		head -n1 age_mask.xy >> age_mask.xy
		;;
	W)	echo "Masking out to the WEST of the trench"
		echo "${WLONA} ${NLATA}" >> age_mask.xy
		echo "${WLONA} ${SLATA}" >> age_mask.xy
 		head -n1 age_mask.xy >> age_mask.xy
		;;
	esac
	
	gmt grdmask age_mask.xy -I1m/1m -R$ABOUNDS -N1/1/NaN -r -Gage_mask.grd
	gmt grdsample age_scaled.grd -R$ABOUNDS -I1m/1m -r -nl -Gage_tmp2.grd
	gmt grdmath age_mask.grd age_tmp2.grd MUL = age_tmp3.grd
	gmt grdcut age_tmp3.grd -R$BOUNDS -Gage_start.grd
	
	# FILTERING
	
	# extend the grid
# 	gmt grdfilter age_start.grd -D3 -Ni -Fg80 -Gagef1.grd
	
	# reset original ages
# 	gmt grdmath age_start.grd agef1.grd AND = agef1_reset.grd
	
	# filter again
# 	gmt grdfilter agef1_reset.grd -D3 -Ni -Fg40 -Gagef2.grd -V
	
	# reset again
# 	gmt grdmath age_start.grd agef2.grd AND = age_tmp4.grd
	
	# now set the NaNs to 1 Ma
# 	gmt grdmath age_tmp4.grd 1 AND = agef2_reset.grd
	
	# filter again but over 20 km. now we have a complete age grid
# 	gmt grdfilter agef2_reset.grd -D3 -Ni -Fg20 -Gage.grd -V
	
	# CLEAN UP FILES
	# rm traxis.xy age.nc age_scaled.grd age_tmp*.grd