#!/usr/bin/perl

use strict;

                                         #  PROGRAMA  #
###########################################################################################################
###########################################################################################################
#### Amb aquest programa predim el gen amb màxima puntuació en l'strand positiu tenint en compte que ##### 
#### coincideixin el remainder i el frame dels exons consecutius, que els exons no se solapin i que   #####
#### el gen compleixi el model de gen i les distàncies entre exons que l'usuari entrarà mitjançant un ####
#### fitxer de text. El programa consta bàsicament  de 4 parts: 

########################    1- ENTRADA DE DADES                             ###############################
########################    2- ALGORITME                                    ###############################
########################    3- SORTIDA DE DADES(en format gff)              ###############################
########################    4- FUNCIONS (que s'usen al llarg del programa)  ############################## 
###########################################################################################################
###########################################################################################################




######################################## 1-ENTRADA DE DADES ###############################################
###########################################################################################################

######Declarem el fitxer GFF sobre el que treballarà el nostre programa  

 my $fitxer;
 
 $fitxer = $ARGV[0];


### Definim un file handle amb el fitxer que ens entren i crearem la matriu on cada exó serà una fila de 
### la matriu i cada columna del fitxer serà una columna de la matriu. És a dir, cada columna 
### de la matriu serà una de les característiques de l'exó.(acceptor, score, donor, strand...)


 open(DADES,"<$fitxer") || die "No es pot obrir el fitxer\n";


 my $i;
 my @Exons;  ### És la matriu on posarem tots els exons ###
 my $remainder;

 $i=0;



 while(<DADES>) {

   
    if( m/(\w+)\s+([\w\.]+)\s+(\w+)\s+(\d+)\s+(\d+)\s+([\d\.]+)\s+(\+)\s+([012])/) { 
	
	$Exons[$i][0]=$1;   ### Entrada de la seqüència a NCBI o Genbank.Normalment serà NM_00000  ###

        $Exons[$i][1]=$2;   ### Referència del programa amb el qual han estat predits els exons ###

	$Exons[$i][2]=$3;   ### Tipus d'exó. Pot ser first, internal o terminal ###

	$Exons[$i][3]=$4;   ### Posició de l'acceptor de l'exó ###
	
        $Exons[$i][4]=$5;   ### Posició del donor de l'exó ###

	$Exons[$i][5]=$6;   ### Puntuació assignada a l'exó. Pot ser positiva o negativa.

	$Exons[$i][6]=$7;   ### Strand o cadena on està predit l'exó. Ha de ser positiu (així ho especifica 
                            ### l'expressió regular). ###

	$Exons[$i][7]=$8;   ### Frame o pauta de lectura de l'exó ###
 
        $Exons[$i][8]= (3 - ((($Exons[$i][4] - $Exons[$i][3] + 1) - $Exons[$i][7]) %3))%3;  ### Crearem una
                                              ### nova columna on posarem el remainder de cada exó calculat

 	

        $i=$i+1;
	

    }

 } 

 close(DADES);


###Ordenem la matriu @Exons per odre creixent d'acceptor: ###

 @Exons= sort{${$a}[3] <=> ${$b}[3]}(@Exons) ; ### La columna 3 de la matriu és la que fa referència a l'acc                                               ###  eptor.

                             
########################################## 2- ALGORITME ### ###############################################
###########################################################################################################



 my $i;
 my $j;
 my @puns;  ### Vector on guardem les puntuacions màximes dels gens ###
 my @gens;  ### Vector on guardem els exons que formen part del millor gen ###

 $i=0;

 while ($i < scalar(@Exons)) {       

   
 $j = 0;
 $puns[$i] = 0;
 $gens[$i] = -1;
   


    while ($j < $i) {

	  if ( $Exons[$j][4] < $Exons[$i][3] && $Exons[$j][8] == $Exons[$i][7] &&
             restriccio ($Exons[$j][2],$Exons[$i][2],($Exons[$i][3] - $Exons[$j][4]))) {  

 	        if ($puns[$j] > $puns[$i]) {

	       	 $puns[$i] = $puns[$j];
	       	 $gens[$i] = $j;
	    
	         }      

	  }

	  $j = $j+1;

       
    } 
 
    $puns[$i] = $puns[$i]+ $Exons[$i][5];
   
    $i = $i+1;
 
    print "la puntuació és @puns\n" ;

 }  

 print "@gens\n";

 print  "la puntuació màxima és "; print max(@puns),"\n";




######################################## 3- SORTIDA DE DADES ##############################################
###########################################################################################################

### Imprimirem les característiques de tots els exons que formen el millor gen ensamblat per tal de deixar 
### ja preparades les dades pel format tabular GFF.###

 my @entrada; 
 my @ref;
 my @tipus;
 my @acceptors;
 my @donors;
 my @score;
 my @cadena;
 my @pauta;

 my $i;
 my $j;

 $i = posmax(@puns);
 $j = 0;

 while ($i >= 0 ) {
    
    $entrada[$j] = $Exons[$i][0];
    $ref[$j] = $Exons[$i][1]; 
    $tipus[$j] = $Exons[$i][2];
    $acceptors[$j] = $Exons[$i][3];
    $donors [$j] = $Exons[$i][4];
    $score[$j] = $Exons[$i][5];
    $cadena[$j] = $Exons[$i][6];
    $pauta[$j] = $Exons[$i][7];

    $i=$gens[$i];

    $j = $j + 1; 

 }

 @entrada=  reverse(@entrada);
 @ref =  reverse(@ref);
 @tipus =  reverse(@tipus);
 @acceptors =  reverse(@acceptors);
 @donors = reverse(@donors);
 @score =  reverse(@score);
 @cadena =  reverse(@cadena);
 @pauta =  reverse(@pauta);


 my $i;
 $i = 0;
 
 print "\n La predicció de PREGEN és:\n";
 while ($i < scalar(@acceptors)) {
    

 print "$entrada[$i] \t pregen \t $tipus[$i] \t $acceptors[$i] \t $donors[$i] \t $score[$i] \t $cadena[$i] \t $pauta[$i]\n";
    
 $i = $i+1;
 
 }

print "\n La puntuació del gen és "; print max(@puns),"\n";

print "\n El número d'exons que formen el gen és: $i\n";

########################################### 4- FUNCIONS  ##################################################
###########################################################################################################
                       

### funció que ens compara el tipus d'exó!

 sub restriccio {

 ###Declarem el fitxer que crearà cada usuari del programa

 my $param = $ARGV[1];
 my $itipus = 0;
 my @tipus;


        open(PARAMETRES,"<$param") || die "No es pot obrir el fitxer\n";

        while (<PARAMETRES>) {

               if (m/(\w+)\s+(\w+)\s+(\d+)\s+(\d+)/) {
	
	       $tipus[$itipus][0] = $1;
     
               $tipus[$itipus][1] = $2;
     
	       $tipus[$itipus][2] = $3;
     
               $tipus[$itipus][3] = $4; 
	 
               $itipus = $itipus+1;
 
               } 

    
         }
 
         close (PARAMETRES);

         my $tipus5p = $_[0];
         my $tipus3p = $_[1];
         my $distancia = $_[2];    
         my $r = 0;
         my $restriccio = 0;
  

         while ( $r < scalar (@tipus) && $restriccio == 0) {

	        if ( lc($tipus5p) eq $tipus[$r][0] && lc($tipus3p) eq $tipus[$r][1]  && $distancia >= $tipus                   [$r][2] && $distancia <= $tipus[$r][3]) {

    	         $restriccio = 1;

	        }

	        $r = $r + 1;
          }

          return $restriccio;
 
 }


### funció max!

 sub max {
 my $i = 0;
 my $max = -99999;

          while ($i < scalar(@_)) {

                if ($_[$i] >= $max) {
      
                $max = $_[$i];
                }

                $i = $i + 1;

           }

           return $max;

 }



### funció posició màxima!

 sub posmax {
 my $i = 0;
 my $max = -99999;
 my $posmax;

            while ($i < scalar(@_)) {

                  if ($_[$i] >= $max) {
     
                  $max = $_[$i];
                  $posmax = $i;
      
                  }  

                  $i = $i + 1;

             }
  
             return $posmax;

 }