# Now we are back to the main working directory
workDir=$1
queryFile=$2
# Make directory for the uniprot taxonomy results
if [ ! -d "$workDir/uniprot_taxonomy_using_orfs" ]; then
 mkdir -p  $workDir/uniprot_taxonomy_using_orfs
fi

# Perform blastx of the predicted ORFs against the virus/phage uniprot database reference (TREMBL)
# Use ORF fasta file (from glimmer3)
blastp -query $queryFile -out $workDir/uniprot_taxonomy_using_orfs/blastx_trembl_glimmer_total_orfs.txt -db /blastdb/uniprot_virus_and_phage_TrEMBL/uniprot_virus_and_phage_TrEMBL_db -outfmt 6 -num_threads 40 -max_target_seqs 1 -evalue 1e-6

# Move blastx output files to a specific directory
if [ ! -d "$workDir/uniprot_taxonomy_using_orfs/blastx_raw_results" ]; then
 mkdir -p $workDir/uniprot_taxonomy_using_orfs/blastx_raw_results
fi

mv $workDir/uniprot_taxonomy_using_orfs/blastx_trembl_glimmer_total_orfs.txt $workDir/uniprot_taxonomy_using_orfs/blastx_raw_results

# Get the gene and contig IDs from the blastx results (output 6; tab delimited file)
# Do this all in a loop to make life easier
if [ ! -d "mkdir $workDir/uniprot_taxonomy_using_orfs/virome_phage_blastx_formatted_for_taxa" ]; then
 mkdir -p $workDir/uniprot_taxonomy_using_orfs/virome_phage_blastx_formatted_for_taxa
fi
for file in $(ls $workDir/uniprot_taxonomy_using_orfs/blastx_raw_results); do
    # Set file name variable
    NAME=$(echo ${file})
    #Assign contig IDs to tmp list file
    cut -f 1 $workDir/uniprot_taxonomy_using_orfs/blastx_raw_results/${file} | sed 's/_.*//' > $workDir/uniprot_taxonomy_using_orfs/virome_phage_blastx_formatted_for_taxa/tmp1_${file}
    #Assign gene hit IDs to another tmp list
    cut -f 2 $workDir/uniprot_taxonomy_using_orfs/blastx_raw_results/${file} | sed 's/.*|.*|\(.*\)/\1/' > $workDir/uniprot_taxonomy_using_orfs/virome_phage_blastx_formatted_for_taxa/tmp2_${file}
    #Paste these files together for reference
    paste $workDir/uniprot_taxonomy_using_orfs/virome_phage_blastx_formatted_for_taxa/tmp1_${file} $workDir/uniprot_taxonomy_using_orfs/virome_phage_blastx_formatted_for_taxa/tmp2_${file} > $workDir/uniprot_taxonomy_using_orfs/virome_phage_blastx_formatted_for_taxa/tmp_paste_${file}
    #Add in the taxonomy data from the uniprot reference database
    echo Assigning taxonomy to ${file}...
    if [[ $NAME =~ .*trembl.* ]]; then
        awk 'FNR==NR { a[$1]=$2"\t"$3"\t"$4"\t"$5"\t"$6; next } $2 in a { print $1"\t"a[$2] }' /PhageHostPredction/HumanSkin/references/UniProt-Virus-Phage/uniprot_reference_phage_and_virus_TrEMBL_taxonomy_table.txt $workDir/uniprot_taxonomy_using_orfs/virome_phage_blastx_formatted_for_taxa/tmp_paste_${file} > $workDir/uniprot_taxonomy_using_orfs/virome_phage_blastx_formatted_for_taxa/tmp_awk_cat_${file}
        echo $NAME contains trembl!
    elif [[ $NAME =~ .*swissprot.* ]]; then
        echo $NAME contains swissprot!
    else
        echo $NAME does not contain trembl or swissprot!
    fi
    #Add @ symbol delimiter at first tab for easier perl parsing
    sed 's/\t/@/' $workDir/uniprot_taxonomy_using_orfs/virome_phage_blastx_formatted_for_taxa/tmp_awk_cat_${file} > $workDir/uniprot_taxonomy_using_orfs/virome_phage_blastx_formatted_for_taxa/format_for_perl_${file}
    #Filter out the contigs that did not have at least one orf match per 10kb.
    #Return the numbers of orfs that had assigned taxonomy to each contig.
    #sed 's/@.*/_/' ./uniprot_taxonomy_using_orfs/virome_phage_blastx_formatted_for_taxa/format_for_perl_${file} | sort | uniq -c | sed 's/^ *//' | sed 's/ /\t/' > ./uniprot_taxonomy_using_orfs/virome_phage_blastx_formatted_for_taxa/tmp_orf_count_${file}
    #Assemble a list of contig numbers, only including those contigs which had at least 1 orf per 10kb
    #awk 'FNR==NR { a[$1]=$2; next } $2 in a { print $2"\t"$1"\t"a[$2]"\t"10000*$1/a[$2] }' ./contig_stats/contig_length_without_greater_sign.txt ./uniprot_taxonomy_using_orfs/virome_phage_blastx_formatted_for_taxa/tmp_orf_count_${file} | awk '$4 > 1' > ./uniprot_taxonomy_using_orfs/virome_phage_blastx_formatted_for_taxa/tmp_awk_list_less_10kb_${file}
    #Get only the contig ID numbers, clean away the other information from each row
    #cut -f 1 ./uniprot_taxonomy_using_orfs/virome_phage_blastx_formatted_for_taxa/tmp_awk_list_less_10kb_${file} > ./uniprot_taxonomy_using_orfs/virome_phage_blastx_formatted_for_taxa/tmp_awk_list_less_10kb_clean_${file}
    #Remove the undersores from the ends of the contig IDs so that grep can match the entire words
    #sed 's/\t/_\t/' ./uniprot_taxonomy_using_orfs/virome_phage_blastx_formatted_for_taxa/tmp_awk_cat_${file} > ./uniprot_taxonomy_using_orfs/virome_phage_blastx_formatted_for_taxa/tmp_awk_cat_no_underscore_${file}
    #Use grep to get the contig rows that match the list of contigs to keep
    #grep -w --file=./uniprot_taxonomy_using_orfs/virome_phage_blastx_formatted_for_taxa/tmp_awk_list_less_10kb_clean_${file} ./uniprot_taxonomy_using_orfs/virome_phage_blastx_formatted_for_taxa/tmp_awk_cat_no_underscore_${file} > ./uniprot_taxonomy_using_orfs/virome_phage_blastx_formatted_for_taxa/coverage_filtered_${file}
    #Again add @ symbol delimiter at first tab for easier perl parsing on filtered contigs
    #sed 's/_\t/@/' ./uniprot_taxonomy_using_orfs/virome_phage_blastx_formatted_for_taxa/coverage_filtered_${file} > ./uniprot_taxonomy_using_orfs/virome_phage_blastx_formatted_for_taxa/format_for_perl_filtered_${file}
    # Remove the strain specific IDs at the ends of the phage and virus species names (for exmple, staph_phage_0594 is just staph_phage)
    # Here I am primarily changing the phage identifiers and am not changing much of the other taxonomic information
    # There were some duplicates that we noticed upon manual inspection, so I also fix those here
    sed 's/\(\t.*_phage\).*$/\1/' $workDir/uniprot_taxonomy_using_orfs/virome_phage_blastx_formatted_for_taxa/format_for_perl_${file} | sed 's/\(\t.*_virus\).*$/\1/' | sed 's/_type.*$//' | sed 's/\(virus\)[^\t^\n]*$/\1/' | sed 's/_bacteriophage.*$/_phage/' | sed 's/\(_phage\).*_phage.*$/\1/' | sed 's/\([A-Za-z]\{2\}phage\)_.*$/\1/' | sed 's/_prophage.*$/_phage/' | sed 's/-/_/' | sed 's/\./_/' | sed 's/Streptomyce_phage/Streptomyces_phage/' | sed 's/Mycobacteriophage/Mycobacterium_phage/' | sed 's/Environmental_Halophage/Environmental_halophage/' | sed 's/Corynephage/Corynebacterium_phage/' | sed 's/Enterobacterial_phage/Enterobacteria_phage/' | sed 's/Enterobacterio_phage/Enterobacteria_phage/' > $workDir/uniprot_taxonomy_using_orfs/virome_phage_blastx_formatted_for_taxa/trimmed_format_for_perl_filtered_${file}
    # Use the uniprot orf contig taxonomy script for contig taxonomy assignment
    perl /PhageHostPredction/Prophage/Virsorter/contig_id_by_orfs.pl $workDir/uniprot_taxonomy_using_orfs/virome_phage_blastx_formatted_for_taxa/trimmed_format_for_perl_filtered_${file} $workDir/uniprot_taxonomy_using_orfs/virome_phage_blastx_formatted_for_taxa/perl_taxonomy_results_${file}
    # Add underscore to the end of each contig number in the file
    sed 's/\t/_\t/' $workDir/uniprot_taxonomy_using_orfs/virome_phage_blastx_formatted_for_taxa/perl_taxonomy_results_${file} > $workDir/uniprot_taxonomy_using_orfs/virome_phage_blastx_formatted_for_taxa/perl_taxonomy_results_underscore_${file}
done