# Giant Virus Finder, a software package to find giant viruses in metagenomes.
# Copyright (C) 2015, Csaba Kerepesi

# This file is part of Giant Virus Finder.

# Giant Virus Finder 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
# 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: <http://www.gnu.org/licenses/>.

# For any other inquiries send an Email to Csaba Kerepesi: kerepesi@caesar.elte.hu
# Giant Virus Finder website: http://pitgroup.org/giant-virus-finder

# $ARGV[0]: target dir
# $ARGV[1]: input metagenome
# $ARGV[2]: cutoff
# Run: perl GiantVirusFinder.pl Giants-in-16-soil-metagenome mgm4477803.3.050.upload.fastq.fasta 0.2

$cutoff=1-$ARGV[2];
open (IN, "TopGiants-300000-allacc.csv");
while (<IN>) {
    chomp;
    my @split_line=split(/\|/,$_);
    $is_giant_acc{$split_line[0]}=1;
}
close IN;
open IN, "$ARGV[0]/phase1-results/Candidates-in-$ARGV[1]";
while (<IN>) {
    if ($_ =~ /^>(.+)/) {
        $act_title=$1;
        $read{$1}=$_;
    }else {
        $read{$act_title}="$read{$act_title}$_";
    }
}
close IN;
open (IN, "$ARGV[0]/phase2-results/reblast-$ARGV[1].blast");
while (<IN>) {
    if ($_ =~ /^Query= (.+)/ ) {
        $act_query=$1;
        $is_topsc_query=1;
        $hit{$act_query}=$read{$1};
    }elsif ($_ =~ /^>(\S+)/) {
        $act_subject_acc=$1;
        $hit{$act_query}="$hit{$act_query}$_";
        $topsc_subject=1;
    }elsif ($_ =~ / Score =\s+(\S+)/) {
        if ($is_topsc_query) {
            $act_top_score=$1;
            if (!exists $is_giant_acc{$act_subject_acc}) {
                $is_giant_query{$act_query}=0;
            }else {
                $is_giant_query{$act_query}=1;
                $topsc_query{$act_query}=$act_top_score;
            }
            $is_topsc_query=0;
        }elsif (($topsc_subject) && ($1>=$cutoff*$act_top_score) && (!exists $is_giant_acc{$act_subject_acc})) {
            $is_giant_query{$act_query}=0;
        }
        if ($topsc_subject) {
            $hit{$act_query}="$hit{$act_query}$_";
        }    
        $topsc_subject=0;
    }elsif ($_ =~ /^Length=/) {
        $hit{$act_query}="$hit{$act_query}$_";
    }
}
close IN;
open OUT2,">","$ARGV[0]/phase2-results/reblast-$ARGV[1].blast-GiantVirusFinder-$ARGV[2].fasta";
open OUT,'>',"$ARGV[0]/phase2-results/reblast-$ARGV[1].blast-GiantVirusFinder-$ARGV[2]-with_hits.txt";
foreach my $key (sort { $topsc_query{$b} <=> $topsc_query{$a} } keys %topsc_query) {
    if ($is_giant_query{$key}) {
        print OUT "$hit{$key}\n";
        print OUT2 $read{$key};
    }
}
close OUT;
close OUT2;
print "$ARGV[0]/phase2-results/reblast-$ARGV[1].blast-GiantVirusFinder-$ARGV[2]-with_hits.txt\n";
print "$ARGV[0]/phase2-results/reblast-$ARGV[1].blast-GiantVirusFinder-$ARGV[2].fasta\n";
