#!/usr/bin/perl -w

use strict;



#CercaGen pretèn ser un programa d'enamblatge d'exons.
#Requereix com a entrada un fitxer amb format .gff amb tots els exons que es vulguin avaluar i un fitxer amb el model de gen.
#Busca la concatenació d'exons de manera que la puntuació resultant sigui màxima.
#Els exons ecadenats són compatibles en la pauta de lectura, no es superposen i compleixen el model de gen introduït.

#COS DEL PROGRAMA

#En cas que el nombre d'arguments passats al programa sigui incorrecte dona error i acaba.

if(scalar(@ARGV)!=3){
    print "\nNOMBRE DE PARAMETRES INCORRECTE: CercaGen.pl Entrada.gff ModelGen.txt Sortida.gff\n\n";
    exit(0);
}

my @Epos=llegir_fitxer($ARGV[0],'+');           #Llegeix els exons de cadena positiva del fitxer d'entrada i els guarda en una matriu
my @Eneg=llegir_fitxer($ARGV[0],'-');           #LLegeix els exons de cadena negativa del fitxer d'entrada i els guarda en una matriu
my @M=llegir_fitxer_model($ARGV[1]);            #LLegeix el fitxer del model de gen a utilitzar i ho guarda en una matriu.
@Epos=ordena(\@Epos,'a');                       #Ordena la matriu d'exons de cadena positiva per odre creixent de acceptor.
@Eneg=ordena(\@Eneg,'a');                       #Ordena la matriu d'exons de cadena negativa per odre creixent de acceptor.
my @Rpos=BestGene(\@Epos,\@M);                  #Calcula el millor gen amb exons de cadena positiva  i torna el resultat en una matriu
my @Mneg=crea_model_neg(\@M);                   #Crea el model de gen que serà utilitzat per la cadena negativa.
my @Rneg=BestGene(\@Eneg,\@Mneg);               #Calcula el millor gen amb exons de cadena negativa i torna el resultat en una matriu

mostra_matriu(\@Rpos);                          #Mostra per pantalla la matriu de resultats per a la cadena positiva.
print "\n";
mostra_matriu(\@Rneg);                          #Mostra per pantalla la matriu de resultats per a  la cadena negativa.
imprimeix_fitxer(\@Rpos,\@Rneg,$ARGV[2],\@Epos,\@Eneg);



#FUNCIO:LLEGIR_FITXER:Llegeix un fitxer .gff i posa les dades en una matriu.
#ARGUMENTS:El nom del fitxer, i la cadena a llegir.
#RETORNA:La matriu amb les dades.

sub llegir_fitxer{
    
    my $i=0;
    my @M;

    if(open(FILE,"<$_[0]")){                     
                                                #Mentre es llegeix el fitxer va guardant totes les dades en una matriu
	while(<FILE>) {                         #Expressió regular que captura camps separats per tabuladors (format .gff)
	                                         
	    if (m/(\A[^#\t]+)\t([^\t]+)\t([^\t]+)\t([^\t]+)\t([^\t]+)\t([^\t]+)\t([^\t]+)\t([^\t\n]+)\n*\t*([^\t\n]*)\n*/){   

		      if($7 eq $_[1]){          #Només agafa els exons amb la cadena corresponent al segon parametre de la funció
		  
		      $M[$i][0]=$1;
		      $M[$i][1]=$2;
		      $M[$i][2]=$3;
		      $M[$i][3]=$4;
		      $M[$i][4]=$5;
		      $M[$i][5]=$6;
		      $M[$i][6]=$7;
		      $M[$i][7]=$8;
		      $M[$i][8]=$9
;
		      $i++;}
	     }
	}
	close(FILE);

	return @M;                             #Torna una matriu plena amb les dades del fitxer.
    }
    else{
                                               #Si no es pot obrir el fitxer dona error i acaba el programa
	print "NO S'HA POGUT OBRIR EL FITXER";
	exit(0);
    }
}


#FUNCIO:LLEGIR_FITXER_MODEL:Llegeix el fitxer del model de gen i guarda les dades en una matriu.
#PARAMETRES:El nom del fitxer.
#RETORNA: La matriu amb les dades.

sub llegir_fitxer_model{

    my $i=0;
    my @M;

    if(open(FILE,"<$_[0]")){

	while(<FILE>){                        #Capturem les dades del fitxer mitjançant una expressió regular.

	    if (m/(\w+)\s+(\w+)\s+(\d+)\s+(\d+)/){
		
		$M[$i][0]=lc($1);             #Guardem el nom del tipus d'exò en lletra minúscula.(funció lc)
		$M[$i][1]=lc($2);
		$M[$i][2]=$3;
		$M[$i][3]=$4;
		$i++;
	    }
	}
	close(FILE);
    }
                                              #Si no es pot obrir el fitxer dona error i acaba el programa
    else{
                      
	print "NO S'HA POGUT OBRIR EL FIXER";
	exit(0);
    }
    return @M;
}


#FUNCIÓ:MOSTRA_MATRIU:Imprimeix per pantalla la matriu de dades
#ARGUMENTS:La matriu de dades per referència
#RETORNA:Res.

sub mostra_matriu{

    my $i;
    my $j;
    my @M=@{$_[0]};

    for($i=0;$i<scalar(@M);$i++){
	print "$i ";                         #Per a cada exò imprimeix també a quina posició de la matriu es troba.
	for($j=0;$j<9;$j++){
	    print"$M[$i][$j]\t";
	}
	print reminder($i,\@M);              #Per a cada exó imprimeix també el seu reminder.
	print "\n";
    }
}

#FUNCIO:ORDENA:Ordena una matriu segons els "donnors" o "acceptors" dels exons que contè.Fa una ordenació per SELECCIÓ
#ARGUMENTS:La matriu de dades passada per referència, "d" o "a" segons com s'hagi d'ordenar.
#RETORNA:La matriu ordenada.

sub ordena{

    my @M=@{$_[0]};
    my $x;
    my $i;
    my $j;
    my $index_menor;

    if($_[1]eq'a'){                          #Ordenem segons acceptor
	$x=3;
    }
    else{
	$x=4;                                #Ordenem segons donnor
    }

    for($i=0;$i<scalar(@M)-1;$i++){          

	$index_menor=$i;

	for($j=$i+1;$j<scalar(@M);$j++){

	    if($M[$j][$x]<$M[$index_menor][$x]){

		$index_menor=$j;
	    }
	}
	@M=intercanvia($i,$index_menor,\@M);
    }
    return @M;
}    


#FUNCIO:INTERCANVIA:Intercanvia dues files d'una matriu.És una funció necessària per a la funció ORDENA.
#ARGUMENTS:Les dues posicions a intercanviar i la matriu.
#RETORNA:La matriu resultant de l'operació.

sub intercanvia{

    my @temp;
    my $i;
    my @M=@{$_[2]};
    
    for($i=0;$i<9;$i++){
	$temp[$i]=$M[$_[0]][$i];
	$M[$_[0]][$i]=$M[$_[1]][$i];
	$M[$_[1]][$i]=$temp[$i];
    }
    return @M;
}
     

#FUNCIO:BESTGENE:Calcula el millor gen a partir de la matriu d'exons.
#PARAMETRES:La matriu amb les dades dels exons ordenada per acceptor i la matriu amb el model de gen.
#RETORNA:Una matriu amb les dades dels exons que conformen el millor gen possible

sub BestGene{

    my @E=@{$_[0]};           
    my $i;
    my $j;
    my @S;                            #Vector amb la puntuació màxima per a cada gen acabat en cada exó.
    my @G;                            #Vector amb els gens que ha usta per construir cada gen.
    my @R;                            #Matriu de resultats.
    

    for($i=0;$i<scalar(@E);$i++){
	
	$S[$i]=0;
	$G[$i]=-1;

	for($j=0;$j<=$i-1;$j++){  

#L'exò incorporat no s'ha de solapar amb l'anterior, ha de millorar la puntuació global,ha d'entrar en fase (frame igual al reminder de     l'exò anterior) i ha de complir les restriccions imposades pel model de gen.
	   
	    if(($E[$j][4]<$E[$i][3])&&($S[$j]>$S[$i])&&(reminder($j,\@E)==$E[$i][7])&&model_gene($_[1],$E[$j][2],$E[$i][2],$E[$i][3]-$E[$j][4])){
	       		
		$S[$i]=$S[$j];
		$G[$i]=$j;
	
	    }
	    
	}
	$S[$i]=$S[$i]+$E[$i][5];
    }

    return mostra_gen(\@G,\@S,\@E);   #Crida la funció que generarà la matriu complerta de resultats.
    #mostra_matriu(\@R);
    
}

#FUNCIO:MAXIM_VECTOR:Busca el valor màxim d'un vector.
#PARAMETRES:El vector.
#RETORNA: La posició màxima.

sub maxim_vector{

    my @V=@{$_[0]};
    my $i;
    my $cont=0;
    my $pos;

    for($i=0;$i<scalar(@V);$i++){

	if($V[$i]>$cont){

	    $cont=$V[$i];
	    $pos=$i;
	}
    }
    return $pos;
}
    
    
#FUNCIÓ:MOSTRA_GEN: Crea un vector amb els exons que composen el gen amb major puntuació.
#PARAMETRES:Els vectors de puntuacions i gens usats creats per BestGenr
#RETORNA:Un vector amb les posicions dels exons.
   
sub mostra_gen{
    
    my @G=@{$_[0]};
    my @S=@{$_[1]};
    my @E=@{$_[2]};
    my $i=0;
    my @R;

    my $x=maxim_vector(\@S);
 
    
    while($x!=-1){

 	#print "$E[$x][3]-\>$E[$x][4]_____$x\n";
	$R[$i]=$x;
	$x=$G[$x];
	$i++;
    }
    @R=reverse(@R);                   
   
    return crea_resultat(\@R,\@E);            #Crida una funció que crearà la matriu de dades resulat.
}

#FUNCIO:CREA_RESULTAT:Crea una matriu amb només els exons que formen el millor gen.
#ARGUMENTS:La matriu de dades i el vector amb les posicions del exons del millor gen.
#RETORNA: La matriu de resultats

sub crea_resultat{

    my @E=@{$_[1]};
    my @R=@{$_[0]};
    my @M;
    my $i;
    my $j;
    my $k=0;
    my $x;

    while ($k<scalar(@R)){
	$x=$R[$k];
	for($i=0;$i<9;$i++){
	   
	    $M[$k][$i]=$E[$x][$i];
		
	    }
	
    $k++;
    }


return @M;
}


    
#FUNCIO:REMINDER:Calcula el reminder d'un exó.
#ARGUMENTS:L'exó.
#RETORNA:El reminder.

sub reminder{

    my @E=@{$_[1]};
    my $j=$_[0];

    return (3-(($E[$j][4]+1-$E[$j][3]-$E[$j][7])%3))%3;
 
}    


#FUNCIO:IMPRIMEIX_FITXER:Crea el fitxer de resultat, imprimeix la matriu de resultats i afegeix cometraris. Té format .gff
#PARAMETRES:Les matrius de resultats i de dades per a les dues cadenes(positiva i negativa)
#RETORNA:Res.

sub imprimeix_fitxer{

    my @Rpos=@{$_[0]};
    my @Rneg=@{$_[1]};
    my @Epos=@{$_[3]};
    my @Eneg=@{$_[4]};
    my $i;
    my $j;
    my $data=localtime(time);                        #La data en que es crea el fitxer
    my $scorepos=0;
    my $scoreneg=0;
    my $longpos=0;
    my $longneg=0;

    for($i=0;$i<scalar(@Rpos);$i++){                 #Calcula la longitud i puntuació sumada de tots els exons  
	$scorepos+=$Rpos[$i][5];
	$longpos+=$Rpos[$i][4]-$Rpos[$i][3]+1;
    }
    for($i=0;$i<scalar(@Rneg);$i++){
	$scoreneg+=$Rneg[$i][5];
	$longneg+=$Rneg[$i][4]-$Rneg[$i][3]+1;
    }
		
    
    open(RES,"> $_[2]");
    
    print RES "#CercaGen v 1.0\n";
    print RES "#$data\n";
    print RES "#Sequence $Epos[0][0].\n";
    print RES "#Found Forward Exons: ",scalar(@Rpos),"\t\tScore: $scorepos\tLenght: $longpos bps\n";
    print RES "#Found Reverse Exons: ",scalar(@Rneg),"\t\tScore: $scoreneg\tLenght: $longneg bps\n";
    
    print RES "\n";
    
    for($i=0;$i<scalar(@Rpos);$i++){

	for($j=0;$j<9;$j++){
	    print RES "$Rpos[$i][$j]\t";
	}
	print RES "\n";
    }
    print RES "\n";
    for($i=0;$i<scalar(@Rneg);$i++){

	for($j=0;$j<9;$j++){
	    print RES "$Rneg[$i][$j]\t";
	}
	print RES "\n";
    }
    
    close(RES);
}


#FUNCIO:MODEL_GENE:Comprova si els exons que se'nsamblaràn es corresponen amb el model de gen.
#ARGUMENTS: El tipus dels dos exons, i la distància que hi ha entre ells.
#RETORNA: 0 si no són compatibles i 1 si ho són.

sub model_gene {

    my $i;
    my @ModelGen=@{$_[0]};
    my $tipus1=$_[1];
    my $tipus2=$_[2];
    my $dist=$_[3];
    my $res=0;


    for($i=0;$i<scalar@ModelGen&&!$res;$i++){
	
	if(lc($tipus1)eq $ModelGen[$i][0] && lc($tipus2)eq $ModelGen[$i][1] && $dist>=$ModelGen[$i][2] && $dist<$ModelGen[$i][3]){

	    $res=1;
	}
    }
    return $res;
}



#FUNCIO:CRA_MODEL_NEG:Adapta la matriu de model de gen a l'strand negatiu.
#ARGUMENTS:La matriu amb el model de Gen.
#RETORNA: Una nova matriu adaptada a la cadena negativa.

sub crea_model_neg{

    my @Mneg;
    my @Mpos=@{$_[0]};
    my $i;
    my $temp;
    
    for($i=0;$i<scalar(@Mpos);$i++){

	$Mneg[$i][0]=$Mpos[$i][1];
	$Mneg[$i][1]=$Mpos[$i][0];
	$Mneg[$i][2]=$Mpos[$i][2];
	$Mneg[$i][3]=$Mpos[$i][3];
    }

    return @Mneg;
}


#Aquest programa té dues limitacions important:

                                               #Funciona sota un ordre de complexitat quadràtic.
                                               #Tracta completament per separat la cadena negativa i la positiva.

