#!/bin/perl
use strict;
use warnings;
use Data::Dumper;

#Variablen und Pfade der Ursprungsdaten und CoFoldDaten initialisieren

my $UrsprungsDaten_ct="/u/home/praktikum/RNA_SandraGerstl/RNA_STRAND_data/LongRNAs/";
my $CoFoldDaten_ct="/u/home/praktikum/RNA_SandraGerstl/RNA_STRAND_data/LongRNAs/CofoldStructure/bpseq_to_ct/";

#Verzeichnisse/Listen mit Daten angelegt, nur ct-Files sollen verwendet werden
#Hash wird mit key (Dateinname ohne Endung) und Pfad der Dateien gefuellt
   
#Tabellenkopf
print "Key\tSequenzlaenge\tSensitivity\tPPV\tMCC\tSpecificity\n";


my %Ursprungsdaten_pfade=();
opendir(my $UrsprungsDaten,$UrsprungsDaten_ct) || die "Kann Verzeichnis Referenzdaten nicht oeffnen!";
my @UrsprungsDaten_Liste=grep {/\.ct$/} readdir($UrsprungsDaten);
closedir $UrsprungsDaten;
for my $eintrag (@UrsprungsDaten_Liste) {
    my $key= $eintrag;
    $key=~s/\.ct$//i;
    $Ursprungsdaten_pfade{$key}="$UrsprungsDaten_ct"."$eintrag";
}

my %CoFoldDaten_pfade=();
opendir(my $CoFoldDaten,$CoFoldDaten_ct) || die "Kann Verzeichnis CoFoldDaten nicht oeffnen!";
my @CoFoldDaten_Liste=grep {/\.ct$/} readdir($CoFoldDaten);
closedir $CoFoldDaten;
for my $eintrag2 (@CoFoldDaten_Liste) {
    my $key2= $eintrag2;
    $key2=~s/\.ct$//i;
    $CoFoldDaten_pfade{$key2}="$CoFoldDaten_ct"."$eintrag2";


}
#Vergleich der hashs %Ursprungsdaten_pfade und %CoFoldDaten_Pfade


for my $key ( keys %CoFoldDaten_pfade) {
    my $CoFoldDaten_pfad=$CoFoldDaten_pfade{$key};
    my $Ursprungsdaten_pfad=$Ursprungsdaten_pfade{$key};

    my $TP=0;
    my $FP=0;
    my $TN=0;
    my $FN=0;

    #print "$key\n";
    #print "$CoFoldDaten_pfad\n";
    #print "$Ursprungsdaten_pfad\n\n";

    open(IN,"<$CoFoldDaten_pfad");
    my @CoFoldDaten_Vergleichsliste = (<IN>);

    open(IN,"<$Ursprungsdaten_pfad");
    my @Ursprungsdaten_Vergleichsliste = (<IN>);

    #print "$Ursprungsdaten_pfad\n";
    #print "$CoFoldDaten_pfad\n";

    my $seq_laenge=0;
    my $line_diff= grep {/^#/} @Ursprungsdaten_Vergleichsliste;
    for (my $i=1; $i<=$#CoFoldDaten_Vergleichsliste;$i++) {
	$seq_laenge++;
	my @zeile_cofold=split(/\s+/, $CoFoldDaten_Vergleichsliste[$i]);
	
	my @zeile_ursprungsdaten = split (/\s+/, $Ursprungsdaten_Vergleichsliste[$i+$line_diff]);

        #test nach korrekten Datensaetzen mit entsprechenden Zeilen

	if ($zeile_cofold[1] != $zeile_ursprungsdaten[1]){
	    die "Zeilen entsprechen nicht den korrekten Zeilen im zugeordneten File!";
	}

	#count TP
	if ($zeile_cofold[5] == $zeile_ursprungsdaten[5] && $zeile_cofold[5] > $zeile_cofold[1]) {
	    $TP++;
	}

	
	#count FP
		if ($zeile_cofold[5] != $zeile_ursprungsdaten[5] && $zeile_cofold[1] < $zeile_cofold[5]) {
	    $FP++;
	}
	
	#count FN
		if ($zeile_cofold[5] != $zeile_ursprungsdaten[5]  && $zeile_ursprungsdaten[1] < $zeile_ursprungsdaten[5]) {
	    $FN++;
	}

    }
   
    #Total berechnen
    my $total;
    $total = (($seq_laenge - 3)*($seq_laenge - 4))/2;
    
    #TN berechnen
    $TN= $total - $TP - $FP - $FN;
    
    #print stats
    #print "Datei: $key, True Positives: $TP, True Negatives:$TN, False Negatives: $FN, False Positives: $FP, Total=$total\n";

    #Sensitivity
    my $sensitivity = $TP / ($TP + $FN);
    my$ger_sensitivity = sprintf("%.6f", $sensitivity);

    #PPV
    my $PPV = $TP / ($TP + $FP);
    my$ger_PPV = sprintf("%.6f", $PPV);

    #MCC
    my $MCC = (($TP * $TN) - ($FP * $FN)) /  ((($TP + $FP) * ($TP + $FN) * ($TN + $FP) * ($TN + $FN)) **0.5 );
    my$ger_MCC = sprintf("%.6f", $MCC);

    #Specificity
    my $specificity = $TN / ($FP + $TN);
    my$ger_specificity = sprintf("%.6f", $specificity);


    #print "Sensitivity: $sensitivity, PPV: $PPV, MCC: $MCC \n";
    print "$key\t$seq_laenge\t$ger_sensitivity\t$ger_PPV\t$ger_MCC\t$ger_specificity\n";
    #In Tabelle einspeichern
    


}