#!/usr/bin/perl -w

use strict;



#########################################################################################################################################
#########################################################################################################################################
#########################################################################################################################################
###                                                                                                                                   ###
###                                                                                                                                   ###
###                                                                                                                                   ###
###                                                         SARCORFS 2004                                                             ###
###                                                                                                                                   ###
###                                                                                                                                   ###
###                                                                                                                                   ###
###                                                                                                                                   ###
###                                                                                                                                   ###
###                                         Ana Cuadrado Garcia & David Dominguez Fandos                                              ###
###                                                                                                                                   ###
###                                                                                                                                   ###
###                                                                                                                                   ###
###                                            Facultat de Ciencies de la Salut i la Vida                                             ###
###                                                                                                                                   ###
###                                                     Universitat Pompeu Fabra                                                      ###
###                                                                                                                                   ###
###                                                                                                                                   ###
###                                                                                                                                   ###
#########################################################################################################################################
#########################################################################################################################################
#########################################################################################################################################




##############################################DECLARACION DE LAS VARIABLES A UTILIZAR#####################################################

my @entrada;                              # vector en el cual en cada posicion hay un nucleotido de $entrada_seguida
my @complementaria;                       # vector con la cadena complementaria, cada posicion es un codon de la cadena complementaria
my $j;                                    # variable para recorrer el vector @complementaria para generar la secuencia complementaria seguida
my $entrada_seguida;                      # variable que contiene la secuencia problema introducida y su complementaria sin saltos de linea. 
my $pauta;                                # variable que especifica la pauta de lectura en la cual el programa lee $entrada_seguida
my $direccion;                            # variable que especifica en que direccion el programa lee la secuencia (forward o reverse)
my @entrada_codones;                      # vector en el cual registramos en cada posicion uno de los codones de $entrada_seguida
my $c;                                    # variable que recorre el vector @entrada, para generar @complementaria
my $identificador;                        # variable en la cual introducimos la primera linea de la secuencia problema
my @lineas_entrada;                       # vector que contiene en cada posicion una linea de la secuencia problema, excepto la primera que se encuentra en $identificador
my $s;                                    # variable que recorre las posiciones del vector @lineas_entrada
my $p;                                    # variable que recorre las posiciones de @entrada_codones
my $q;                                    # variable que recorre las posiciones de @orfs
my @orfs;                                 # vector que guarda en cada posicion un codon del ORF encontrado (que empiece por ATG y tenga 60 o mas aminoacidos)
my %codon_aminoacido;                     # hash con los diferentes codones del codigo genetico y los aminoacidos correspondientes a cada codon
my $m;                                    # variable que recorre el vector @orfs cuando el ORF tiene una longitud mayor o igual a 60 aminoacidos
my $clave;                                # variable que contiene cada vez un codon diferente de @orfs (palabra clave del hash) y asi poder obtener el valor que se asocia a $clave                                       
my $l;                                    # variable que recorre cada posicion del vector @complementaria para construir la cadena complementaria de la secuencia problema


#creamos un hash con los codones del codigo genetico (%codon_aminoacido) en el cual la palabra clave es el codon y el valor asociado es el aminoacido correspondiente en el codigo de una letra. El simbolo "*" significa que el codon es un codon stop.

  %codon_aminoacido= ("TTT" => "F",
                      "TTC" => "F",
                      "TTA" => "L",
                      "TTG" => "L",
                      "CTT" => "L",
                      "CTC" => "L",
                      "CTA" => "L",
                      "CTG" => "L",
                      "ATT" => "I",
                      "ATC" => "I",
                      "ATA" => "I",
                      "ATG" => "M",
                      "GTT" => "V",
                      "GTC" => "V",
                      "GTA" => "V",
                      "GTG" => "V",
                      "TCT" => "S",
                      "TCC" => "S",
                      "TCA" => "S",
                      "TCG" => "S",
                      "CCT" => "P",
                      "CCC" => "P",
                      "CCA" => "P",
                      "CCG" => "P",
                      "ACT" => "T",
                      "ACC" => "T",
                      "ACA" => "T",
                      "ACG" => "T",
		      "GCT" => "A",
		      "GCC" => "A",
		      "GCA" => "A",
		      "GCG" => "A",
		      "TAT" => "Y",
		      "TAC" => "Y",
		      "TAA" => "*",
		      "TAG" => "*",		    
		      "CAT" => "H",
		      "CAC" => "H",
		      "CAA" => "Q",
		      "CAG" => "Q",
		      "AAT" => "N",
		      "AAC" => "N",
		      "AAA" => "K",
		      "AAG" => "K",
		      "GAT" => "D",
		      "GAC" => "D",
		      "GAA" => "E",
		      "GAG" => "E",
		      "TGT" => "C",
		      "TGC" => "C",
		      "TGA" => "*",
		      "TGG" => "W",
		      "CGT" => "R",
		      "CGC" => "R",
		      "CGA" => "R",
		      "CGG" => "R",
		      "AGT" => "S",
		      "AGC" => "S",
		      "AGA" => "R",
		      "AGG" => "R",
		      "GGT" => "G",
		      "GGC" => "G",
		      "GGA" => "G",
		      "GGG" => "G");





# construimos una composicion alternativa (o secuencial) que nos pida que le introduzcamos el fichero con la secuencia problema en formato FASTA en caso de no haberla introducido.



if (scalar(@ARGV) < 1)  {

    print "Falta introducir el fichero con la secuencia a analizar en formato FASTA\n";

    exit;

}



open (RESULTADOS, ">orfs_conseguidos.fa");  # redireccionamos a un fichero en formato FASTA los resultados obtenidos con la ejecucion del programa

funcion_entrada_seguida($ARGV[0]);  # llamamos a la funcion funcion_entrada_seguida


########################################################################################################################################
#######################################################PAUTAS DE LECTURA FORWARD#########################################################
#########################################################################################################################################



##############################################Direccion forward, Pauta de lectura 0######################################################


$direccion = "Forward";

$pauta = 0;

@entrada_codones = ($entrada_seguida =~ m/[ACGT][ACGT][ACGT]/g);   # encajamos la expresion regular en $entrada_seguida tantas veces como se pueda y lo vamos introduciendo en el vector @entrada_codones

traduccion_orfs();                                                 # llamamos a la funcion traduccion_orfs



##############################################Direccion forward, Pauta de lectura 1#####################################################


$entrada_seguida = substr ($entrada_seguida,1);        # sustraemos de la secuencia recogida en $entrada_seguida un nucleótido para poder leer en la pauta de lectura 1

$direccion = "Forward";

$pauta = 1;

@entrada_codones = ($entrada_seguida =~ m/[ACGT][ACGT][ACGT]/g);

traduccion_orfs();


##############################################Direccion forward, Pauta de lectura 2#####################################################


$entrada_seguida = substr ($entrada_seguida,1);        # sustraemos un nuevo nucleotido para obtener la pauta de lectura 2

$direccion = "Forward";

$pauta = 2;

@entrada_codones = ($entrada_seguida =~ m/[ACGT][ACGT][ACGT]/g);

traduccion_orfs();


#########################################################################################################################################
##############################################OBTENCION DE LA CADENA COMPLEMENTARIA######################################################
#########################################################################################################################################


funcion_complementaria();    # llamamos a la funcion funcion_complementaria

funcion_complementaria_seguida();  # llamamos a la funcion funcion_complementaria_seguida


#########################################################################################################################################
####################################################PAUTAS DE LECTURA REVERSE############################################################
#########################################################################################################################################



##############################################Direccion reverse, Pauta de lectura 0######################################################


$direccion = "Reverse";

$pauta = 0;

@entrada_codones = ($entrada_seguida =~ m/\w\w\w/g);

traduccion_orfs();

##############################################Direccion reverse, Pauta de lectura 1######################################################


$entrada_seguida = substr ($entrada_seguida,1);        # sustraemos de la secuencia recogida en $entrada_seguida un nucleotido para poder leer en la pauta de lectura 1

$direccion = "Reverse";

$pauta = 1;

@entrada_codones = ($entrada_seguida =~ m/\w\w\w/g);

traduccion_orfs();


##############################################Direccion reverse, Pauta de lectura 2######################################################


$entrada_seguida = substr ($entrada_seguida,1);       # sustraemos un nuevo nucleotido para obtener la pauta de lectura 2

$direccion = "Reverse";

$pauta = 2;

@entrada_codones = ($entrada_seguida =~ m/\w\w\w/g);

traduccion_orfs();





#########################################################################################################################################
#############################################ESPECIFICACION DE FUNCIONES UTILIZADAS######################################################
#########################################################################################################################################



############################################OBTENCION DE LA CADENA FORWARD SEGUIDA#######################################################



sub funcion_entrada_seguida {

###NOMBRE: funcion_entrada_seguida
###PARAMETRO: $ARGV [0], argumento que contiene la secuencia problema FASTA a analizar
###PROPOSITO: obtener la secuencia introducida en formato FASTA toda seguida (sin saltos de linea), e introducirla en la variable $entrada_seguida
###RESULTADO QUE DEVUELVE: la variable $entrada_seguida; variable que contendra toda la secuencia inicial en una sola linea. Sin la linea del identificador.

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

$identificador = <FILE>;         # le pasamos a la variable $identificador la primera linea de la secuencia problema a analizar, que contiene el identificador

@lineas_entrada = <FILE>;        # vector que contiene en cada posicion una linea de la secuencia problema, excepto la primera que se encuentra en $identificador


# a continuacion construimos una composicion iterativa para recorrer todas las posiciones del vector @lineas_entrada, el contenido de las cuales ira a parar a $entrada_seguida para construir la secuencia problema seguida 

$s = 0;                          # variable que recorre las posiciones del vector @lineas_entrada

$entrada_seguida = "";           # variable que contendra toda la secuencia problema seguida (sin saltos de linea)

while ($s < scalar(@lineas_entrada))   {

    chomp ($lineas_entrada[$s]);

    $entrada_seguida = $entrada_seguida . $lineas_entrada[$s];   # concatenacion de las diferentes lineas de la secuencia FASTA en $entrada_seguida

    $s = $s + 1;

}


close (FILE);                     # cerramos el fichero ya que se deja de utilizar

return $entrada_seguida;          # nos devuelve el resultado que ha calculado la funcion en este caso $entrada_seguida

}


##############################################SINTESIS DE LA CADENA COMPLEMENTARIA#######################################################
 

sub funcion_complementaria  {

###NOMBRE: funcion_complementaria
###PROPOSITO: construir la cadena complementaria de la secuencia introducida en formato FASTA, a partir de la secuencia registrada en $entrada_seguida
###RESULTADO QUE DEVUELVE: la cadena complementaria introducida en el vector @complementaria, en el cual cada posicion es un nucleotido de la secuencia


@entrada = split (//,$entrada_seguida);  #funcion que permite cortar la secuencia almacenada en $entrada_seguida y la introduce en el vector @complementaria donde encontraremos un simbolo por posicion

$c= scalar(@entrada)-1;                  #empieza por la ultima posicion, por lo tanto la longitud del vector menos uno

$l = 0;

while ($c >= 0) {                        ######bucle para hacer la complementaria de la cadena almacenada en $entrada_seguida#####

    if ($entrada[$c] eq 'C') {    
	$complementaria[$l] = 'G';
    }

    if ($entrada[$c] eq 'G') {
	$complementaria[$l] = 'C';
    }

    if ($entrada[$c] eq 'A') {
	$complementaria[$l] = 'T';
    }

    if ($entrada[$c] eq 'T') {
       	$complementaria[$l] = 'A';
    }

    $c = $c - 1;                         # variable para que continue en la siguiente posicion de @entrada
    $l = $l + 1;                         # variable para avanzar en el vector @complementaria
} 

return @complementaria;

}     


########################################OBTENCION DE LA CADENA COMPLEMENTARIA SEGUIDA####################################################



sub funcion_complementaria_seguida   {

###NOMBRE: funcion_complementaria_seguida
###PROPOSITO: obtener la cadena complementaria de la secuencia introducida en formato FASTA sin saltos de linea, e introducirla en la variable $entrada_seguida
###RESULTADO QUE DEVUELVE: la variable $entrada_seguida; variable que contendra toda la secuencia complementaria en una sola linea

$j = 0;                      # variable que recorre el vector @complementaria

$entrada_seguida = "";       # inicialmente $entrada_seguida esta vacia, no contiene ningun elemento

# construimos una composicion iterativa para recorrer todas las posiciones de @complementaria, el contenido de las cuales ira a parar a $entrada_seguida para construir la secuencia complementaria seguida

while ($j < scalar (@complementaria)) {

    chomp ($complementaria [$j]);

    $entrada_seguida = $entrada_seguida . $complementaria[$j];   # concatenaciones parciales para construir $entrada_seguida (en este caso la complementaria)

    $j = $j + 1;

}

return $entrada_seguida;

}


###########################BUSQUEDA DE METIONINAS Y TRADUCCION DE OPEN READING FRAMES MEDIANTE UNA FUNCION###############################

sub traduccion_orfs   {


###NOMBRE: traduccion_orfs
###PROPOSITO: encontrar las metioninas (ATG) en @entrada_codones y traducir los ORFs si tienen 60 o mas aminoacidos
###RESULTADO QUE DEVUELVE: los ORFs traducidos en el fichero orfs_conseguidos.fa (formato FASTA)


$p = 0;                                    # variable que recorre las posiciones de @entrada_codones
                                    

# construimos una composicion iterativa para recorrer todas las posiciones del vector @entrada_codones y de esta manera que encuentre ATG (principio de ORF)


  while ($p < scalar (@entrada_codones))   {  

      $q = 0;                             # variable que recorre las posiciones de @orfs
   
    if ($entrada_codones[$p] eq "ATG")   {

# construimos otra composicion iterativa que guarda los ORF en el vector @orfs mientras no encuentre ningun codon stop

	while (($p < scalar (@entrada_codones)) && ($entrada_codones[$p] ne "TAG") && ($entrada_codones[$p] ne "TAA") && ($entrada_codones[$p] ne "TGA"))  {

	    $orfs[$q] = $entrada_codones[$p];
	    
	    $p = $p+1;
	    $q = $q+1;
	}                                                     # llave que cierra el while que tiene arriba ($p < scalar (@entrada_codones)
    
        if ($q >= 60)   {

            print RESULTADOS "> $direccion, Pauta:$pauta, Inicio-Final del ORF:",($p*3-$q*3)+1+$pauta,"-",($p*3)+$pauta,"\n";

            $m = 0;                                           # variable que recorre el vector @orfs cuando el ORF tiene una longitud >= 60 aminoacidos

            while ($m < $q)   {

               $clave = $orfs[$m];

               print RESULTADOS "$codon_aminoacido{$clave}";   # se obtienen las claves del hash correspondientes a cada codón

               $m = $m + 1;
   
            }                                                  # llave que cierra el while ($m < $q)

            print RESULTADOS "\n\n\n";

       }                                                       # llave que cierra el if ($q >= 60)
      
    }                                                          # llave que cierra el primer if ($entrada_codones[$p] eq "ATG")
   
   $p = $p + 1;

  }                                                            # llave que cierra el primer while ($p < scalar (@entrada_codones))

}                                                              # llave que cierra la funcion


#########################################################################################################################################

print "Se ha generado el fichero orfs_conseguidos.fa donde se han almacenado los resultados obtenidos con la ejecucion del programa SARCORFS\n";





#########################################################################################################################################
####################################FINAL DEL PROGRAMA QUE BUSCA Y TRADUCE ORFS SI TIENE >= 60 AMINOACIDOS###############################
#########################################################################################################################################







