hlapers

#!/usr/bin/env bash

# script directory
if [[ "$OSTYPE" == "linux-gnu" ]]; then

    myexecpath="$( readlink -f "$0" )"

elif [[ "$OSTYPE" == "darwin"* ]]; then
        
    myexecpath="$( readlink "$0" )"

else
    myexecpath="$( pwd )"
fi

execdir="$( dirname "$myexecpath" )"
scriptdir=$execdir/script

# Usage functions
usage() {
    echo "Usage: hlapers [modes]"
    echo ""
    printf "%-20s %s\n" "prepare-ref" "Prepare transcript fasta files."
    printf "%-20s %s\n" "index" "Create index for read alignment."
    printf "%-20s %s\n" "bam2fq" "Convert BAM to fastq."
    printf "%-20s %s\n" "genotype" "Infer HLA genotypes."
    printf "%-20s %s\n" "quant" "Quantify HLA expression."
}

usage_prepare() {
    echo "Usage: hlapers prepare-ref [options]"
    echo ""
    printf "%-20s %s\n" "-t | --transcripts" "Fasta with Gencode transcript sequences."
    printf "%-20s %s\n" "-a | --annotations" "GTF from Gencode for the same Genome version."
    printf "%-20s %s\n" "-i | --imgt" "Path to IMGT directory."
    printf "%-20s %s\n" "-o | --out" "Output directory."
}

usage_index() {
    echo "Usage: hlapers index [options]"
    echo ""
    printf "%-20s %s\n" "-t | --transcripts" "Fasta with Gencode transcript sequences."
    printf "%-20s %s\n" "-p | --threads" "Number of threads."
    printf "%-20s %s\n" "-o | --out" "Output directory."
    printf "%-20s %s\n" "--kallisto" "Create index for kallisto pipeline instead of STARsalmon."
}

usage_genotype() {
    echo "Usage: hlapers genotype [options]"
    echo ""
    printf "%-20s %s\n" "-i | --index" "Index generated by 'hlapers index'."
    printf "%-20s %s\n" "-t | --transcripts" "Fasta with Gencode transcripts sequences used for 'hlapers index'." 
    printf "%-20s %s\n" "-1 | --fq1" "Fastq for READ 1."
    printf "%-20s %s\n" "-2 | --fq2" "Fastq for READ 2."
    printf "%-20s %s\n" "-p | --threads" "Number of threads."
    printf "%-20s %s\n" "-o | --outprefix" "Output prefix name."
    printf "%-20s %s\n" "--kallisto" "Use kallisto for genotyping."
}

usage_quant() {
    echo "Usage: hlapers quant [options]"
    echo ""
    printf "%-20s %s\n" "-t | --transcripts" "Reference transcripts directory."
    printf "%-20s %s\n" "-g | --genotypes" "*_genotypes.tsv file generated by 'hlapers genotype'."
    printf "%-20s %s\n" "-1 | --fq1" "Fastq for READ 1."
    printf "%-20s %s\n" "-2 | --fq2" "Fastq for READ 2."
    printf "%-20s %s\n" "-p | --threads" "Number of threads."
    printf "%-20s %s\n" "-o | --out" "Output prefix name."
    printf "%-20s %s\n" "--salmonreads" "Use Salmon lightweight alignment for quantification (NOT TESTED)"
    printf "%-20s %s\n" "--kallisto" "Use kallisto for quantification."
}

usage_bam2fq() {
    echo "Usage: hlapers bam2fq [options]"
    echo ""
    printf "%-20s %s\n" "-m | --mhc-coords" "Genomic coordinates of the MHC region in chrN:start-end format if MHC fastq is desired."
    printf "%-20s %s\n" "-b | --bam" "BAM file (if -m is specified, needs to be sorted by coordinate; otherwise use --sort)."
    printf "%-20s %s\n" "-o | --outprefix" "Output prefix name."
    printf "%-20s %s\n" "--sort" "Sort input BAM file by coordinate (REQUIRED if -m is specified and BAM is not sorted by coordinate)."
}

# Mode functions
run_prepref () {
    
    if [[ "$#" -eq 0 ]]; then
	usage_prepare
	exit 1
    fi
    
    while [[ "$#" -gt 0 ]]; do
        case "$1" in
        -h|--help)
            usage_prepare
            exit
            ;;
        -t|--transcripts)
            transcripts="$2"
            shift 2
            ;;
        -a|--annotations)
            annotations="$2"
            shift 2
            ;;
	-i|--imgt)
	    imgt="$2"
	    shift 2
	    ;;
	-o|-out)
	    out="$2"
	    shift 2
	    ;;
        *)
            echo "ERROR: unknown parameter $1"
            usage_prepare
            exit 1
            ;;
        esac
    done
        
    Rscript $scriptdir/make_index_files.R $transcripts $annotations $imgt $out
}


run_index() {
    
    if [[ "$#" -eq 0 ]]; then
	usage_index
	exit 1
    fi
    
    while [[ "$#" -gt 0 ]]; do
        case "$1" in
        -h|--help)
            usage_index
            exit
            ;;
	-t|--transcripts)
	    transcripts="$2"
	    shift 2
	    ;;
	-p|--threads)
	    threads="$2"
	    shift 2
	    ;;
	-o|--out)
	    out="$2"
	    shift 2
	    ;;
	--kallisto)
	    k=1
	    shift
	    ;;
        *)
            echo "ERROR: unknown parameter $1"
            usage_index
            exit 1
            ;;
        esac
    done
  
    kallisto=${k-0}

    if [[ "$kallisto" -eq 0 ]]; then

	$scriptdir/run_starindex.sh $transcripts $out $threads

    elif [[ "$kallisto" -eq 1 ]]; then
        
        $scriptdir/run_kallistoindex.sh $transcripts $out

    fi
}

run_genotype() {

    if [[ "$#" -eq 0 ]]; then
        usage_genotype
        exit 1
    fi

    while [[ "$#" -gt 0 ]]; do
	case "$1" in
	    -h|--help)
		usage_genotype
		exit
		;;
	    -i|--index)
		index="$2"
		shift 2
		;;
	    -t|--transcripts)
		transcripts="$2"
		shift 2
		;;
	    -1|--fq1)
		fq1="$2"
		shift 2
		;;
	    -2|--fq2)
		fq2="$2"
		shift  2
		;;
	    -p|--threads)
		threads="$2"
		shift 2
		;;
	    -o|--out)
		outprefix="$2"
		shift 2
		;;
	    --kallisto)
		k=1
		shift
		;;
	    *)
		echo "ERROR: unknown parameter $1"
		usage_genotype
		exit 1
		;;
	esac
    done
    
    kallisto=${k-0}

    if [[ "$kallisto" -eq 0 ]]; then

	$scriptdir/run_genotyping.sh $index $transcripts $fq1 $fq2 $outprefix $threads

    elif [[ "$kallisto" -eq 1 ]]; then

	$scriptdir/run_genotypingkallisto.sh $index $fq1 $fq2 $outprefix $threads

    fi

}

run_quant() {
    if [[ "$#" -eq 0 ]]; then
        usage_quant
        exit 1
    fi

    while [[ "$#" -gt 0 ]]; do
	case "$1" in
	    -h|--help)
		usage_quant
		exit
		;;
	    -t|--transcripts)
		transcripts="$2"
		shift 2
		;;
	    -g|--genotypes)
		genos="$2"
		shift 2
		;;
	    -1|--fq1)
		fq1="$2"
		shift 2
		;;
	    -2|--fq2)
		fq2="$2"
		shift 2
		;;
	    -p|--threads)
		threads="$2"
		shift 2
		;;
	    -o|--outprefix)
		outprefix="$2"
		shift 2
		;;
	    --kallisto)
		k=1
		shift
		;;
	    --salmonreads)
		s=1
		shift
		;;
	    *)
		echo "ERROR: unknown parameter $1"
		usage_quant
		exit 1
		;;
	esac
    done
    
    kallisto=${k-0}
    salmonreads=${s-0}

    if [[ "$kallisto" -eq 0 ]] && [[ "$salmonreads" -eq 0 ]]; then

	$scriptdir/run_quant.sh $transcripts $genos $fq1 $fq2 $outprefix $threads
    
    elif [[ "$kallisto" -eq 1 ]] && [[ "$salmonreads" -eq 0 ]]; then

	$scriptdir/run_quantkallisto.sh $transcripts $genos $fq1 $fq2 $outprefix $threads

    elif [[ "$kallisto" -eq 0 ]] && [[ "$salmonreads" -eq 1 ]]; then

	$scriptdir/run_quantSalmonReads.sh $transcripts $genos $fq1 $fq2 $outprefix $threads

    fi
}

run_bam2fq() {
    if [[ "$#" -eq 0 ]]; then
        usage_bam2fq
        exit 1
    fi

    while [[ "$#" -gt 0 ]]; do
	case "$1" in
	    -h|--help)
		usage_bam2fq
		exit
		;;
	    -b|--bam)
		bam="$2"
		shift 2
		;;
	    -m|--mhc)
		if [[ -n "$2" ]] && [[ "$2" != -* ]]; then
		    mhc="$2"
		    shift 2
		else
		    echo "ERROR: missing MHC coordinates."
		    exit 1
		fi
		;;
	    -o|--out)
		out="$2"
		shift 2
		;;
	    --sort)
		s=1
		shift
		;;
	    *)
		echo "ERROR: unknown parameter $1"
		usage_bam2fq
		exit 1
		;;
	esac
    done

    sort=${s-0}

    if [[ -n "$mhc" ]]; then
	$scriptdir/run_bam2fq_mhc.sh $bam $mhc $out $sort
    else
	$scriptdir/run_bam2fq.sh $bam $out $sort
    fi

}

# main
if [[ "$#" -eq 0 ]]; then
    usage
    exit 1
fi

while [[ $# -gt 0 ]]
do
    case "$1" in
    -h|--help)
	usage
	exit
	;;
    
    prepare-ref)
	shift
	run_prepref "$@"
	break
	;;

    index)
	shift
	run_index "$@"
	break
	;;

    genotype)
	shift
	run_genotype "$@"
	break
	;;

    quant)
	shift
	run_quant "$@"
	break
	;;
    bam2fq)
	shift
	run_bam2fq "$@"
	break
	;;
    *)
	echo "ERROR: unknown parameter $1"
	usage
	exit 1
	;;
    esac
done