LABEL

#!/usr/bin/env bash
set -Eeuo pipefail

err_report() {
  echo -e "\e[1m\e[31mError on line $(caller)\e[0m" >&2
}

trap "err_report" ERR INT QUIT ABRT PIPE TERM

# LABEL - Lineage Assignment By Extended Learning
# Predicts the clade for a given nucleotide sequence & gene.

# Affiliation: Centers for Disease Control & Prevention
# Author: Samuel S. Shepard et al. (vfn4@cdc.gov)

PROGRAM="LABEL"
YEAR="2016"
VERSION="v0.5.2"
AUTHOR="Samuel S. Shepard"
EMAIL="vfn4@cdc.gov"
AFFIL="Centers for Disease Control & Prevention"

# GPL version 3
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program.  If not, see <http://www.gnu.org/licenses/>.
#
#
# DICLAIMER & LIMITATION OF LIABILITY.
#
# The materials embodied in this software are "as-is" and without warranty
# of any kind, express, implied or otherwise, including without limitation,
# any warranty of fitness for a particular purpose. In no event shall the
# Centers for Disease Control and Prevention (CDC) or the United States (U.S.)
# Government be liable to you or anyone else for any direct, special, incidental,
# indirect or consequential damages of any kind, or any damages whatsoever,
# including without limitation, loss of profit, loss of use, savings or revenue,
# or the claims of third parties, whether or not CDC or the U.S. Government has
# been advised of the possibility of such loss, however caused and on any theory
# of liability, arising out of or in connection with the possession, use or
# performance of this software.  In no event shall any other party who modifies
# and/or conveys the program as permitted according to GPL license 
# [www.gnu.org/licenses/], make CDC or the U.S. government liable for damages,
# including any general, special, incidental or consequential damages arising out
# of the use or inability to use the program, including but not limited to loss
# of data or data being rendered inaccurate or losses sustained by third parties
# or a failure of the program to operate with any other programs.  Any views,
# prepared by individuals as part of their official duties as United States
# government employees or as contractors of the United States government and
# expressed herein, do not necessarily represent the views of the United States
# government. Such individuals’ participation in any part of the associated work
# is not meant to serve as an official endorsement of the software. The CDC and
# the U.S. government shall not be held liable for damages resulting from any
# statements arising from use of or promotion of the software that may conflict
# with any official position of the United States government.
#


################
# INSTALLATION #
################
# Set manual base path or allow script to locate its current folder automatically.
# The resource folder will be expected to be in the same folder.
# Idea courtesy: stackoverflow.com/questions/59895/can-a-bash-script-tell-what-directory-its-stored-in
bpath=
if [ "$bpath" == "" ]; then
	bpath=$( cd -P "$( dirname "${BASH_SOURCE[0]}" )" && pwd )
fi
resources=LABEL_RES

# OS and ARCH
if [ -x "`which uname`" ]; then
	OS=`uname -s`
	if [ "$(uname -m)" != "x86_64" -a "$OS" == "Linux" ]; then
		OS="Linux32"
	fi
else
	OS="Linux"
fi

## SHELL VARS ##
LANG=POSIX
TAB=$'\t'
shopt -u nocaseglob
shopt -u nocasematch
################
################


# FNC - HELP #
function print_usage() {
	echo -e "\n$PROGRAM $VERSION, updated $YEAR\n$AUTHOR ($EMAIL), $AFFIL"
        echo -e "GPL version 3. This program comes with ABSOLUTELY NO WARRANTY. This is free software."
	echo -e "You are welcome to redistribute it under certain conditions. See:  <http://www.gnu.org/licenses/>.\n\n"
	echo -e "Usage:\n\t$(basename $PROGRAM) [-P MAX_PROC] [-E C_OPT] [-W WRK_PATH|-O OUT_PATH] [-G|-TACRD|-S] [-L LIN_PATH] <nts.fasta> <project> <Module:H5,H9,etc.>"
	echo -e "\t\t-T\tDo TRAINING again instead of using classifier files."
	echo -e "\t\t-A\tDo ALIGNMENT of re-annotated fasta file (sorted by clade) & build its ML tree."
	echo -e "\t\t-C\tDo CONTROL alignment & ML tree construction."
	echo -e "\t\t-E\tSGE clustering option. Use 1 or 2 for SGE with array jobs, else local."
	echo -e "\t\t-R\tNo RECURSIVE prediction. Limits scope, useful with -L option."
	echo -e "\t\t-D\tNo DELETION of extra intermediary files."
	echo -e "\t\t-S\tShow available protein modules."
	echo -e "\t\t-W\tWeb-server mode: requires ABSOLUTE path to WRITABLE working directory."
	echo -e "\t\t-O\tOutput directory path, do not use with web mode."
	echo -e "\t\t-G\tCreate a scoring matrix using given header annotations for Graphing."
	echo -e "Example: $0 -C gisaid_H5N1.fa Bird_Flu H5\n"
}

# FNC - ERROR TEST #
# Idea courtesy: steve-parker.org/sh/exitcodes.shtml
function err_test() {
	if [ "$1" -ne "0" ]; then
		echo ""
		if [ $# -eq "1" ]; then
			echo "`basename $0` ERROR: operations have been ABORTED!"
		else
			echo "`basename $0` ERROR ($2): operations have been ABORTED!"
		fi
		if [ $NO_DELETE -eq 0 -a -d $ppath -a $DO_WEB -eq 0 ]; then
			rm -rf $ppath
		fi
		exit 1
	fi
}

# FNC - CHECK PROGRAM #
# Idea courtesy: stackoverflow.com/questions/592620/check-if-a-program-exists-from-a-bash-script
function check_prgm() {
	command -v $1 &> /dev/null \
	|| { echo "`basename $0` ERROR: Program '$1' not found, please check your PATH or install it." >&2; exit 1; }
}

# PROCESS BASH OPTIONS #
LOCAL_PROC=16	# MAX CONCURRENT PARALLEL
GRID_PROC=128	# MAX CONCURRENT GRID PROCESSORS
DO_TRAIN=0	# Use original data for training.
DO_CONTROL=0	# Control tree computation (for debugging)
NO_RECURSION=0	# Do not use recursion.
NO_DELETE=0	# Do not clean up extra data.
DO_ALIGNMENT=0	# Alignment of query sequences.
DO_WEB=0	# Web-server mode
DO_GRAPH=0
DO_OUT=0
LINEAGE_PATH=""
OUTPUT_DIR=`pwd`
USE_SGE=0	# use SGE loosely
CLUSTER_OPT=0	# clustering option


while getopts "IGP:E:TACRDSW:L:O:h" option;do
	case $option in
		G)
		DO_GRAPH=1
		;;
		P)
		P=$OPTARG
		;;
		E)
		CLUSTER_OPT=$OPTARG
		;;
		T)
		DO_TRAIN=1
		;;
		C)
		DO_CONTROL=1
		;;
		R)
		NO_RECURSION=1
		;;
		D)
		NO_DELETE=1
		;;
		L)
		LINEAGE_PATH=$OPTARG
		;;
		O)
		DO_OUT=1
		OUTPUT_DIR=$OPTARG
		if [ ! -w "$OUTPUT_DIR" ]; then
			echo "`basename $0` ERROR: output directory '$OUTPUT_DIR' not writable."
			exit 1
		fi
		;;
		S)
		ls $bpath/$resources/training_data/
		exit 0
		;;
		A)
		DO_ALIGNMENT=1
		;;
		W)
		DO_WEB=1
		tpath=$OPTARG
		if [ ! -w "$tpath" ]; then
			echo "`basename $0` ERROR: output directory '$tpath' not writable."
			exit 1
		fi
		;;
		h)
		print_usage
		exit 0
		;;
	esac
done

if [ $# -lt $((OPTIND+2)) ]; then
	print_usage
	exit 1
fi

if [ $DO_WEB -eq 1 -a $DO_OUT -eq 1 ]; then
	print_usage
	exit 1
fi


SEQ_LIMIT=100
if [ "$CLUSTER_OPT" -gt 0 ]; then
	USE_SGE=1
fi

# VARIABLES & PATHS #
args=("$@")
inputFasta=${args[$OPTIND-1]}
project=${args[$OPTIND]}
module=${args[$OPTIND+1]}
owd=`pwd`

# Check for GNU coreutils, perl and zip programs
check_prgm zip
check_prgm cut
check_prgm paste
check_prgm perl
check_prgm /usr/bin/env
check_prgm parallel

# Check bundled programs
bin=$bpath/$resources/scripts
hmmscore=$bin/hmmscore_$OS
align2model=$bin/align2model_$OS
modelfromalign=$bin/modelfromalign_$OS
shogun=$bin/shogun_$OS
check_prgm $hmmscore
check_prgm $align2model
check_prgm $modelfromalign
check_prgm $shogun

# Check that 3rd party programs installed
FastTreeMP=FastTreeMP
mafft=mafft
check_prgm $FastTreeMP
check_prgm $mafft


if [ $DO_WEB -eq 0 ]; then
	tpath=$bpath/$resources/test_data
fi

if [ $DO_WEB -eq 1 ]; then
	ppath=$tpath/${project}
else
	ppath=$bpath/$resources/test_data/${project}
fi

if [ -d $ppath ]; then
	echo "`basename $0` ERROR: Project currently in use. See: $ppath"
	err_test 1
else
	mkdir -p $ppath
fi

if [ "$LINEAGE_PATH" != "" ]; then
	tnpath=$bpath/$resources/training_data/$module/$LINEAGE_PATH
	mpath=$tnpath
	grouping=$LINEAGE_PATH

else
	tnpath=$bpath/$resources/training_data/$module
	mpath=$tnpath
	grouping=$module
fi
level=`basename $grouping`

if [ ! -d $tnpath ]; then
	echo -n "`basename $0` ERROR: module '$module' does not exist. Try: LABEL -S"
	err_test 1
fi

# FNC - SCORE_LINEAGE #
# Scores sequences using a particular HMM for a lineage.
function scoreLineage() {
	local ppath=$1
	local m=$2
	local level=$3
	local db=$4
	local run=`basename $m .mod`
	local mpath=`dirname $m`

	cd $ppath
	if [ -r $mpath/null.mod ]; then
		$hmmscore $run -db $db -modelfile $m -subtract_null 0 &> /dev/null
	elif [ "$(basename $ppath)" == "x-rev" ] || [ -d $ppath/x-rev ]; then
		$hmmscore $run -db $db -modelfile $m -dpstyle 1 -subtract_null 1 &> /dev/null
	else
		$hmmscore $run -db $db -modelfile $m &> /dev/null
	fi
	$bin/parseScores.pl $ppath/$run.dist > $ppath/$run.tab
}

# FNC - DO_ALIGNMENT #
# Perform a mafft alignment if available, otherwise use muscle.
function doAlignment() {
	local in=$1
	local out=$2

	if [ $DO_WEB -eq 0 ]; then
		if [ $USE_MAFFT -eq 1 ]; then
			$mafft --maxiterate 1000 --thread -1 $in > $out
		else
			$muscle -in $in -out $out
		fi
	else
		if [ $USE_MAFFT -eq 1 ]; then
			$mafft --maxiterate 1000 --thread -1 $in > $out 2> /dev/null
		else
			$muscle -in $in -out $out &> /dev/null
		fi
	fi
	rm $in
}

# FNC - DO_HMM_ALIGNMENT #
# Performs alignment using SAM's align2model.
function doHMMalignment() {
	local in=$1
	local out=$2
	local mod=$3

	cd $ppath/FASTA
	$align2model $out -db $in -modelfile $mod &> /dev/null
	$bin/removeInsertionStates.pl $out.a2m > $out
	rm $in $out.a2m
}

# FNC - DO_ANALYSIS #
# Main recursive function for prediction at each level.
function doAnalysis() {
	local ppath=$1
	local mpath=$2
	local tnpath=$3
	local grouping=$4
	local level=$5
	local script=$ppath/${level}.tmp.sh
	local m
	local joblist
	local check
	local c

	# Concurrency limiting courtesy "tangens"
	local tmpfile=$ppath/${project}.tab.tmp
	if [ -f $tmpfile ]; then
		rm $tmpfile
	fi
	[ -d $mpath/x-rev -a ! -d $ppath/x-rev ] && mkdir $ppath/x-rev

	local HMMSCORE_OPTS=""
	if [ -r $mpath/null.mod ]; then
		local -a modList=($mpath/null.mod $mpath/*hmm.mod)
		HMMSCORE_OPTS=" -subtract_null 0"
	elif [ -d $mpath/x-rev ]; then
		local -a modList=($mpath/*hmm.mod $mpath/x-rev/*hmm.mod)
		HMMSCORE_OPTS=" -dpstyle 1 -subtract_null 1"
	else
		local -a modList=($mpath/*hmm.mod)
	fi

	# number of modules, sequences, and groups
	local mods=${#modList[@]}
	local seqs=$(grep '>' $ppath/${project}_${level}.fas -c)
	local procs=$LOCAL_PROC
	local doGrid=0
	if [ "$USE_SGE" -eq 1 -a "$seqs" -ge "$SEQ_LIMIT" ]; then
		procs=$GRID_PROC
		doGrid=1
	fi

	local groups=$(expr $seqs / 100 + 1)
	local capacity=$(expr $procs / $mods + 1)
	[ "$groups" -gt "$capacity" ] && groups=$capacity
	local partitions=$(expr $groups \* $mods)

cat >$script  <<EOL
#!/usr/bin/env bash

set -Eeuo pipefail

trap "echo -e '\\e[31m\\e[1m$script error!\\e[0m'" ERR

LANG=$LANG
shopt -u nocaseglob;shopt -u nocasematch

if [ "\$#" -eq "1" ]; then
	ID=\$1
else
	ID=\$SGE_TASK_ID
fi

ID=\$((\$ID - 1))
m=\$((\$ID / $groups))
leaf=\$(printf %04d \$((\$ID % $groups + 1)))
db=$ppath/leaf_\${leaf}.tmpp

mods=(${modList[*]})
mod=\${mods[\$m]}
run="\$(basename \$mod .mod)_\$leaf"

if [[ "\$(basename \$(dirname \$mod))" == "x-rev" ]]; then
	workdir=$ppath/x-rev
else
	workdir=$ppath
fi

cd \$workdir
$hmmscore \$run -db \$db -modelfile \$mod $HMMSCORE_OPTS
exit 0
EOL
	chmod 755 $script
	$bin/interleavedSamples.pl -X tmpp -G $groups $ppath/${project}_${level}.fas $ppath/leaf &> /dev/null
	if [ "$doGrid" -eq "1" ]; then
		# specifying -q all.q may be very important to the success of this call
		qsub $IRMA_QUEUE -t 1-$partitions:1 -N "N${project}-${level}" -sync y -j y -o $script.o $script &> /dev/null
		if [ "$?" -ne "0" ]; then
			(( $DO_WEB )) || echo -e "\n\nWARNING: qsub of '$project-$level' failed. Switching to master node.\n\n"
			parallel --workdir $ppath -j $LOCAL_PROC $script {} ::: $(seq $partitions) &> /dev/null
		else
			rm $script.o
		fi
	else
		parallel --workdir $ppath -j $LOCAL_PROC $script {} ::: $(seq $partitions) &> /dev/null 
	fi
	(( $DO_WEB )) || echo -n '...'

	local mod=""
	local run=""
	local pat=""
	for mod in "${modList[@]}";do
		run=$(basename $mod .mod)
		pat=$(basename $(dirname $mod))

		if [[ "$pat" == "x-rev" ]]; then
			$bin/parseScores.pl $ppath/x-rev/${run}_????.dist > $ppath/x-rev/${run}.tab && rm $ppath/x-rev/${run}_????.dist
		else
			$bin/parseScores.pl $ppath/${run}_????.dist > $ppath/${run}.tab && rm $ppath/${run}_????.dist
		fi
	done
	rm $script $ppath/*.tmpp

	local tmpnullfile=$ppath/null_????.dist
	if [ -f $tmpnullfile ]; then
		rm $tmpnullfile
	fi

	local nllArg=""
	local trnArg=""
	local fltrArg=""
	local vitCRarg=""
	[ -r "$ppath/null.tab" ] 	&& nllArg="-N $ppath/null.tab"
	[ $DO_TRAIN -eq 1 ] 		&& trnArg="-R"
	[ -r "$tnpath/x-filter.txt" ]	&& fltrArg="-D $tnpath/x-filter.txt"
	[ -d $ppath/x-rev ]		&& vitCRarg="-C $ppath/x-rev -F 3"

	$bin/doLABELlevel.pl $ppath $tnpath $trnArg $fltrArg $nllArg $vitCRarg -S $shogun 2> $ppath/${project}.log
	err_test $? $LINENO

	set +o pipefail
	check=`ls -d $mpath/c-* 2>/dev/null | wc -w` || true
	set -o pipefail
	if [ $check -ne 0 -a $NO_RECURSION -eq 0 ]; then
		for c in $mpath/c-*; do
			c=`basename $c`
			set +o pipefail
			check=`cut -f2 $ppath/LEVEL_result.tab | grep "^${c}$" -c` || true
			set -o pipefail
			if [ $check -ne "0" ]; then
				mkdir -p $ppath/$c
				grep "${TAB}${c}$" $ppath/LEVEL_result.tab \
				  | $bin/fastaExtractor.pl $ppath/${project}_${level}.fas -F 1 \
				  > $ppath/$c/${project}_$c.fas

				if [ "$grouping" == "$module" ]; then
					doAnalysis $ppath/$c $mpath/$c $tnpath/$c $c $c &
				else
					doAnalysis $ppath/$c $mpath/$c $tnpath/$c $grouping/$c $c &
				fi
			fi
		done
	fi
	wait

	if [ $NO_DELETE -eq 0 ]; then
		(rm $ppath/*.o &> /dev/null) || true
		rm $ppath/*_hmm.tab $ppath/*dat $ppath/${project}.log
		rm $ppath/${project}_${level}.fas
		if [ -d $ppath/x-rev ]; then
			rm -rf $ppath/x-rev
		fi
	fi
}

# PROCESS INPUT #
(( $DO_WEB )) || echo "LABEL: Processing input data."
if [ -d $inputFasta ]; then
	echo "`basename $0` ERROR: '$inputFasta' was a directory."
	err_test 1
elif [ ! -r $inputFasta ]; then
	echo "`basename $0` ERROR: '$inputFasta' is missing or had unreadable permissions."
	err_test 1
elif [ ! -s $inputFasta -a ! -p $inputFasta ]; then
	echo "`basename $0` ERROR: '$inputFasta' is an empty input file."
	err_test 1
else
	if [ -p $inputFasta ]; then
		IS_PIPE=1
	else
		IS_PIPE=0
	fi

	$bin/stripSequences.pl -F $inputFasta '0-9 :~.-' | $bin/removeByRedundantHeader.pl > $ppath/${project}_$level.fas
	cp $ppath/${project}_$level.fas $ppath/${project}_$level.fas2


	if [ "$USE_SGE" -eq "1" ]; then
		AVG_LEN=$(expr $(grep -v '>' $ppath/${project}_$level.fas|wc -c) / $(grep '>' $ppath/${project}_$level.fas -c) + 1)
		if [ "$AVG_LEN" -le 300 ]; then
			SEQ_LIMIT=1600
		elif [ "$AVG_LEN" -le 1000 ]; then
			SEQ_LIMIT=400
		elif [ "$AVG_LEN" -le 1800 ]; then
			SEQ_LIMIT=100
		else
			SEQ_LIMIT=50
		fi
	fi
fi

# START MAIN #
cd $ppath

# GRAPHING option
if [ $DO_GRAPH -eq 1 ]; then
	modules=`$bin/findByNamePostOrder.pl -C $mpath mod`
	selected=`$bin/reviseTaxa.pl $ppath/${project}_$level.fas -X '_hmm' -M $modules`

	if [ $DO_TRAIN -eq 1 ]; then
		mkdir $ppath/mods
		mpath=$ppath/mods
		$bin/stripSequences.pl -F $owd/$inputFasta '0-9 ' > $mpath/${project}_${level}.tmp
		$bin/partitionTaxa.pl $mpath/${project}_${level}.tmp $mpath
		cd $mpath
		rm $mpath/${project}_${level}.tmp
		for f in $mpath/*; do
			$bin/removeGapColumns.pl $f
			$modelfromalign ${module}_`basename $f .tmp`_hmm -alignfile $f -alphabet DNA 2>/dev/null
		done
		cd $ppath
		modules=`ls $mpath/*mod`
		selected=`ls $mpath/*mod`
		rm $mpath/*tmp $mpath/*weightoutput
	elif [ $DO_CONTROL -eq 1 ]; then
		$bin/reviseTaxa.pl $ppath/${project}_$level.fas -P 'PRED:' > $ppath/${project}_$level.fas2
		for m in $selected; do
			annot=`basename $m _hmm.mod | cut -f2- -d'_'`
			$bin/selectSequences.pl $tnpath/${module}_downsample.fa -I "{$annot}" >> $ppath/${project}_$level.fas2
		done
		mv $ppath/${project}_$level.fas2 $ppath/${project}_$level.fas
	else
		rm $ppath/*fas2
	fi

	echo -n "LABEL: "
	let i=0
	for m in $selected;do
		joblist=($(jobs -p))
		while (( ${#joblist[*]} >= $P ))
		do
			sleep 0.5
			joblist=($(jobs -p))
		done
		scoreLineage $ppath $m $level &
		(( $DO_WEB )) || echo -n "`basename $m _hmm.mod` "
	done
	wait
	echo "DONE."
	$bin/buildDataMatrix.pl $ppath/${project}_analysis.tab $ppath/*.tab
	err_test $? $LINENO

	if [ $NO_DELETE -eq 0 ]; then
		mv $ppath/${project}_analysis.tab $OUTPUT_DIR
	else
		rm $ppath/*dist
		cd $tpath
		zip -q -r ${project}.zip $project
		cd $owd
		mv $tpath/${project}.zip $OUTPUT_DIR
	fi

	rm -rf $ppath
	exit 0
fi


(( $DO_WEB )) || echo -n "LABEL: Performing lineage/clade prediction."
doAnalysis $ppath $mpath $tnpath $grouping $level
cat `$bin/findByNamePostOrder.pl $ppath LEVEL_result.tab` > $ppath/${project}_recursive.tmp
mv $ppath/${project}_$level.fas2 $ppath/${project}_$level.fas
$bin/finalizeResults.pl $ppath/${project}_recursive.tmp $ppath/${project}_${level}.fas $ppath/$project
mv $ppath/${project}_${level}.fas.final $ppath/${project}_${level}.fas
if [ $DO_WEB -eq 0 ]; then
	echo ""
	cat $ppath/${project}_final.txt
fi
cd $owd
(( $DO_WEB )) || cp $ppath/${project}_final.txt $OUTPUT_DIR
cp $bpath/$resources/README.rtf $ppath
(( $DO_WEB )) || $bin/evaluateResults.pl -H -S $ppath/${project}_final.txt

if [ $NO_DELETE -eq 0 ]; then
	rm $ppath/*tmp
fi

# FASTA MANIPULATION #
mkdir -p $ppath/FASTA
mv $ppath/${project}_${level}.fas $ppath/FASTA/${project}_predictions.fas
$bin/reviseTaxa.pl $ppath/FASTA/${project}_predictions.fas -C -D > $ppath/FASTA/${project}_reannotated.fas
$bin/partitionTaxa.pl $ppath/FASTA/${project}_reannotated.fas $ppath/FASTA/ -P ${project}_clade_
err_test $? $LINENO
cat $ppath/FASTA/${project}_clade* > $ppath/FASTA/${project}_reannotated.fas

if [ $DO_ALIGNMENT -eq 1 -a $DO_CONTROL -eq 1 ]; then
	(( $DO_WEB )) || echo "LABEL: Performing alignments."
	if [ -r "$tnpath/x-control/control.mod" ]; then
		doAlignment $ppath/FASTA/${project}_reannotated.fas $ppath/FASTA/${project}_ordered.fasta $module &
		$bin/interleavedSamples.pl -G $LOCAL_PROC $ppath/FASTA/${project}_predictions.fas $ppath/FASTA/${module}_control -X tmp &> /dev/null
		for sample in $ppath/FASTA/${module}_control*tmp;do
			doHMMalignment $sample `basename $sample`.fasta $tnpath/x-control/control.mod &
		done
		wait
		cat $tnpath/x-control/guide.fasta $ppath/FASTA/${module}_control*tmp.fasta > $ppath/FASTA/${module}_control.fasta
		rm $ppath/FASTA/${module}_control*tmp.fasta
	else
		cat $tnpath/${module}_downsample.fa $ppath/FASTA/${project}_predictions.fas > $ppath/FASTA/${module}_control.fa
		doAlignment $ppath/FASTA/${module}_control.fa $ppath/FASTA/${module}_control.fasta $module &
		doAlignment $ppath/FASTA/${project}_reannotated.fas $ppath/FASTA/${project}_ordered.fasta $module &
		wait
	fi

	$bin/sortByAnnotation.pl $ppath/FASTA/${project}_ordered.fasta > $ppath/FASTA/${project}_ordered.fasta2
	mv $ppath/FASTA/${project}_ordered.fasta2 $ppath/FASTA/${project}_ordered.fasta
	if [ $DO_WEB -eq 1 ]; then
		$FastTreeMP -quiet -quote -nt -gtr -gamma -boot 1000 < $ppath/FASTA/${project}_ordered.fasta > $ppath/FASTA/${project}_tree.nwk 2> /dev/null
		$FastTreeMP -quiet -quote -nt -gtr -gamma -boot 1000 < $ppath/FASTA/${module}_control.fasta > $ppath/FASTA/${module}_control.nwk 2> /dev/null
	else
		echo ""
		echo "LABEL: Calculating maximum likelihood trees."
		echo ""
		$FastTreeMP -quote -nt -gtr -gamma -boot 1000 < $ppath/FASTA/${project}_ordered.fasta > $ppath/FASTA/${project}_tree.nwk
		$FastTreeMP -quote -nt -gtr -gamma -boot 1000 < $ppath/FASTA/${module}_control.fasta > $ppath/FASTA/${module}_control.nwk
		echo ""
	fi
else
	# DO NORMAL TREE/ALIGNMENT #
	if [ $DO_ALIGNMENT -eq 1 ]; then
		(( $DO_WEB )) || echo "LABEL: Making ordered alignment & tree"
		doAlignment $ppath/FASTA/${project}_reannotated.fas $ppath/FASTA/${project}_ordered.fasta $module
		$bin/sortByAnnotation.pl $ppath/FASTA/${project}_ordered.fasta > $ppath/FASTA/${project}_ordered.fasta2
		mv $ppath/FASTA/${project}_ordered.fasta2 $ppath/FASTA/${project}_ordered.fasta
		if [ $DO_WEB -eq 1 ]; then
			$FastTreeMP -quiet -quote -nt -gtr -gamma -boot 1000 < $ppath/FASTA/${project}_ordered.fasta > $ppath/FASTA/${project}_tree.nwk 2> /dev/null
		else
			echo ""
			$FastTreeMP -quote -nt -gtr -gamma -boot 1000 < $ppath/FASTA/${project}_ordered.fasta > $ppath/FASTA/${project}_tree.nwk
			echo ""
		fi
	# DO CONTROL TREE/ALIGNMENT #
	elif [ $DO_CONTROL -eq 1 ]; then
		(( $DO_WEB )) || echo "LABEL: Aligning control datasets."

		if [ -r "$tnpath/x-control/control.mod" ]; then
			$bin/interleavedSamples.pl -G $LOCAL_PROC $ppath/FASTA/${project}_predictions.fas $ppath/FASTA/${module}_control -X tmp &> /dev/null
			for sample in $ppath/FASTA/${module}_control*tmp;do
				doHMMalignment $sample `basename $sample`.fasta $tnpath/x-control/control.mod &
			done
			wait
			cat $tnpath/x-control/guide.fasta $ppath/FASTA/${module}_control*tmp.fasta > $ppath/FASTA/${module}_control.fasta
			rm $ppath/FASTA/${module}_control*tmp.fasta
		else
			cat $tnpath/${module}_downsample.fa $ppath/FASTA/${project}_predictions.fas > $ppath/FASTA/${module}_control.fa
			doAlignment $ppath/FASTA/${module}_control.fa $ppath/FASTA/${module}_control.fasta $module
		fi

		if [ $DO_WEB -eq 1 ]; then
			$FastTreeMP -quiet -quote -nt -gtr -gamma -boot 1000 < $ppath/FASTA/${module}_control.fasta > $ppath/FASTA/${module}_control.nwk 2> /dev/null
		else
			echo ""
			echo "LABEL: Calculating maximum likelihood."
			$FastTreeMP -quote -nt -gtr -gamma -boot 1000 < $ppath/FASTA/${module}_control.fasta > $ppath/FASTA/${module}_control.nwk
			echo ""
		fi
	fi
fi

# CLEAN UP #
if [ "$IS_PIPE" -eq "1" ]; then
	cd $tpath
	zip -q -r ${project}.zip $project
	rm -rf $ppath
	cat ${project}.zip
	rm ${project}.zip
elif [ $DO_WEB -eq 1 ]; then
	echo $ppath
else
	cd $tpath
	zip -q -r ${project}.zip $project
	cd $owd
	mv $tpath/${project}.zip $OUTPUT_DIR
	rm -rf $ppath
fi