#!/usr/bin/perl -w
use strict;
print "Hello! \n";
my $count='0';
my $tmcount='0';
my $tp='0';
my $fp='0';
my $tn='0';
my $fn='0';
my $prefirsttm='0';
my $prelasttm='0';
my $coord;
my @firsttmhmm;
my @lasttmhmm;
my $subunit;
my $tmopm='';
my @opms;
my $tmhmmid;

print "What's the name of the file with TMHMM prediction?\n\n";
my $name=<STDIN>;
chomp $name;
print "\n";
while()
 {if (-e $name)
  {
   open FILE, "<$name" or die ("File can't be opened! $!\n\n");
   last;
  }
  else
  {print "File doesn't exist! \n\nPrint another name of your file \n\n";
   $name=<STDIN>;
   chomp $name;
   print "\n";
  }
 }   

while (my $line=<FILE>)
 {if ($line=~m/\slength:\s/i)
   {$count=$';
    chomp $count;
   }
  elsif ($line=~m/\s+tmhelix\s+(\d+)\s+(\d+)/i)
   {
    $prefirsttm=$1;
    $prelasttm=$2;
    push @firsttmhmm, $prefirsttm;
    push @lasttmhmm, $prelasttm;
    $tmcount=$tmcount+$prelasttm-$prefirsttm+1;
   }
  elsif ($line=~m/\#\s+(.+)\s+Number\sof\spredicted/i)
   {$tmhmmid=$1;
   }
 }           

   close FILE, "<$name" or die ("File can't be closed! $!\n\n");

print "What's the name of the OPM file?\n\n";
my $name1=<STDIN>;
chomp $name1;
print "\n";
while()
 {if (-e $name1)
  {
   open FILE, "<$name1" or die ("File can't be opened! $!\n\n");
   last;
  }
  else
  {print "File doesn't exist! \n\nPrint another name of your file \n\n";
   $name1=<STDIN>;
   chomp $name1;
   print "\n";
  }
 }   

while (my $line=<FILE>)
 {if ($line=~m/\*topology\*\s+subunit\s+/i)
   {$subunit=$';
    chomp $subunit;
    $subunit=~s/\s+.+//g;
   }
  else
   {
   }
 }

   close FILE, "<$name1" or die ("File can't be closed! $!\n\n");
   open FILE, "<$name1" or die ("File can't be opened! $!\n\n");

while (my $line=<FILE>)
 {if ($line=~m/\*$subunit\*.+segments:\s+/i)
   {$tmopm=$';
    chomp $tmopm;
   }
  else
   {
   }
 }


   close FILE, "<$name1" or die ("File can't be closed! $!\n\n");

my $prefirstopm;
my $prelastopm;
my @firstopm;
my @lastopm;

@opms=split /,\s/, $tmopm;
foreach my $i (@opms)
 {$prefirstopm=$i;
  $prelastopm=$i;
  $prefirstopm=~s/\d+\(//g;
  $prefirstopm=~s/-\d+\)//g;
  $prelastopm=~s/\d+\(\d+-//g;
  $prelastopm=~s/\)//g;
  push @firstopm, $prefirstopm;
  push @lastopm, $prelastopm;
 }

print "What's the name of the needle alignment file (beetween analysed protein and protein-prototype)?\n\n";
my $name2=<STDIN>;
chomp $name2;
print "\n";
while()
 {if (-e $name2)
  {
   open FILE, "<$name2" or die ("File can't be opened! $!\n\n");
   last;
  }
  else
  {print "File doesn't exist! \n\nPrint another name of your file \n\n";
   $name2=<STDIN>;
   chomp $name2;
   print "\n";
  }
 }   

my $opmid='';

while (my $line=<FILE>)
 {if ($line=~m/\#\s+\d:\s+(.+)\s*/i && $1 eq $tmhmmid)
   {
   }
  elsif ($line=~m/\#\s+\d:\s+(.+)\s*/i && $1 ne $tmhmmid)
   {$opmid=$1;
   }
 }

   close FILE, "<$name2" or die ("File can't be closed! $!\n\n");
   open FILE, "<$name2" or die ("File can't be opened! $!\n\n");

my $tmhmmseq='';
my $opmseq='';
my @seqtmhmm;
my @seqopm;


while (my $line=<FILE>)
 {if ($line=~m/$tmhmmid\s+\d+\s+([\w-]+)\s+\d+/i)
   {$tmhmmseq=$tmhmmseq.$1;
   }
  elsif ($line=~m/$opmid\s+\d+\s+([\w-]+)\s+\d+/i)
   {$opmseq=$opmseq.$1;
   }
 }


my $del="";
@seqtmhmm = split $del, $tmhmmseq;
@seqopm = split $del, $opmseq;

   close FILE, "<$name2" or die ("File can't be closed! $!\n\n");

my $tmopmcount='0';
for (my $t=0; $t<@firstopm; $t++)
 {$tmopmcount=$tmopmcount+$lastopm[$t]-$firstopm[$t]+1;
 }


for (my $l=0; $l<@seqtmhmm; $l++)
 {if ($seqtmhmm[$l] eq '-')
   {for (my $m=0; $m<@firsttmhmm; $m++)
     {if ($firsttmhmm[$m]<$l+1 && $lasttmhmm[$m]>$l+1)
       {$lasttmhmm[$m]=$lasttmhmm[$m]+1;
        $tp=$tp-1;
       }
      elsif ($firsttmhmm[$m]>$l+1)
       {$firsttmhmm[$m]=$firsttmhmm[$m]+1;
        $lasttmhmm[$m]=$lasttmhmm[$m]+1;
       }
     }
   }
  elsif ($seqopm[$l] eq '-')
   {for (my $n=0; $n<@firstopm; $n++)
     {if ($firstopm[$n]<$l+1 && $lastopm[$n]>$l+1)
       {$lastopm[$n]=$lastopm[$n]+1;
        $tp=$tp-1;
       }
      elsif ($firstopm[$n]>$l+1)
       {$firstopm[$n]=$firstopm[$n]+1;
        $lastopm[$n]=$lastopm[$n]+1;
       }
     }
   }
 }



for (my $k=0; $k<@firsttmhmm; $k++)
 {for (my $j=0; $j<@firstopm; $j++)
  {if ($firsttmhmm[$k]<=$firstopm[$j]&&$lasttmhmm[$k]>=$firstopm[$j])
    {if ($lasttmhmm[$k]<=$lastopm[$j])
       {$tp=$tp+$lasttmhmm[$k]-$firstopm[$j]+1;
       }
     elsif ($lasttmhmm[$k]>$lastopm[$j])
       {$tp=$tp+$lastopm[$j]-$firstopm[$j]+1;
       }
    }
   elsif ($firsttmhmm[$k]>$firstopm[$j]&&$firsttmhmm[$k]<=$lastopm[$j])
    {if ($lasttmhmm[$k]<=$lastopm[$j])
       {$tp=$tp+$lasttmhmm[$k]-$firsttmhmm[$k]+1;
       }
     elsif ($lasttmhmm[$k]>$lastopm[$j])
       {$tp=$tp+$lastopm[$j]-$firsttmhmm[$k]+1;
       }
    }
  }
 }                
                                        
$fn=$tmopmcount-$tp; 
$fp=$tmcount-$tp;
$tn=$count-$tmcount-$fn;
my $sens=sprintf("%.2f",$tp/($tp+$fn));
my $spec=sprintf("%.2f",$tn/($tn+$fp));
my $prec=sprintf("%.2f",$tp/($tp+$fp));
my $overpred=sprintf("%.2f",$fp/($fp+$tp));
my $notpred=sprintf("%.2f",$fn/($tn+$fn));


print "The length is $count aa.";
print "\n";
print "The number of predicted transmembrane aa is $tmcount.\n";
print "The number of true positives (TP) is $tp.\n";
print "The number of false positives (FP) is $fp.\n";
print "The number of true negatives (TN) is $tn.\n";
print "The number of false negatives (FN) is $fn.\n";
print "The sensivity is $sens.\n";
print "The specificity is $spec.\n";
print "The precision is $prec.\n";
print "The overprediction is $overpred.\n";
print "The not-enough-prediction is $notpred.\n";

open FILE, ">tm.txt" or die ("File can't be opened! $!\n\n");
print FILE "The length - $count.\n";
print FILE "The number of predicted transmembrane aa - $tmcount.\n";
print FILE "TP - $tp.\n";
print FILE "FP - $fp.\n";
print FILE "TN - $tn.\n";
print FILE "FN - $fn.\n";
print FILE "Sensivity - $sens.\n";
print FILE "Specificity - $spec.\n";
print FILE "Precision - $prec.\n";
print FILE "Overprediction - $overpred.\n";
print FILE "Not-enough-prediction - $notpred.\n";




<STDIN>;