#!/usr/bin/perl
use strict;
#-----------------------
#para concatemers, gera tags e as estatisticas.
my @sizes;
my $minimum_size = 32;
my $maximum_size = 34;

my $input = "";
my $output = "";
my $routput = "";
my $mono = 1;
my $lost = 0;
my $punc = "catg";

my $version = 0.46;
my $vname = "extract_ditags.pl";

if (!$ARGV[0] || $ARGV[0]eq "-h") {
	print "\n\t-min <integer> \t:: set minimum size of the ditags [default:$minimum_size]\n",
	"\t-max <integer> \t:: set maximum size of the ditags [default:$maximum_size]\n",
	"\t-i <FILE> \t:: input file name [default:STDIN]\n",
	"\t-o <FILE> \t:: output file name [default:STDOUT]\n",
	"\t-r <FILE> \t:: report file name [default:STDERR]\n",
	"\t-m [0-2]  \t:: set the recovery of monotags [default:$mono]\n\t\t0 - none\n\t\t1 - evaluate size\n\t\t2 - always \n",
	"\t-l      \t:: activate the recovery of mutated pontuation sites\n",
	"\t-p <sequence> \t:: set the punctuation sequence [default:$punc]\n",
	"\t-h      \t:: display help information\n",
	"\t-v      \t:: display current version\n",
	"\n\t$vname Version $version\n\n";
	exit;
}

if ($ARGV[0] eq "-v") {
	print ("\n\t$vname Version $version\n\n");
	exit;
}
#
##first read command line to collect arguments
#
while (my $args = shift(@ARGV)) {
	#
	#minimun size
	#
	if ($args eq "-min") {
		if ($args = shift(@ARGV)) {
			$minimum_size = $args;
		}
		else {
			print STDERR "MINIMUM SIZE NOT SPECIFIED!\n";
			exit;
		}
	}

	#	
	#maximun size
	#
	elsif ($args eq "-max") {
		if ($args = shift(@ARGV)) {
			$maximum_size = $args;
		}
		else {
			print STDERR "MAXIMUM SIZE NOT SPECIFIED!\n";
			exit;
		}
	}

	#
	#input file
	#
	elsif ($args eq "-i") {
		if ($args = shift(@ARGV)) {
			$input = $args;
		}
		else {
			print STDERR "INPUT NOT SPECIFIED!\n";
			exit;
		}
	}
	#
	#output file
	#
	elsif ($args eq "-o") {
		if ($args = shift(@ARGV)) {
			$output = $args;
		}
		else {
			print STDERR "OUTPUT NOT SPECIFIED!\n";
			exit;
		}
	}
	#
	#output file
	#
	elsif ($args eq "-r") {
		if ($args = shift(@ARGV)) {
			$routput = $args;
		}
		else {
			print STDERR "REPORT OUTPUT NOT SPECIFIED!\n";
			exit;
		}
	}
	#
	#monotags
	#
	if ($args eq "-m") {
		if ($args = shift(@ARGV)) {
			$mono = $args;
		}
		else {
			$mono = 0;
		}
	}
	#
	#lost punctuation
	#
	if ($args eq "-l") {
		$lost = 1;
	}
	#
	#punctuation
	#
	elsif ($args eq "-p") {
		if ($args = shift(@ARGV)) {
			$punc = $args;
		}
		else {
			print STDERR "PUNCTUATUON NOT SPECIFIED!\n";
			exit;
		}
	}

}

if ($input) {
	open INPUT, "<$input" or die "COULD NOT OPEN $input :: $!\n";
	$input = "<INPUT>";
}
else {
	$input = "<>";
}

if ($output) {
	open OUTPUT, ">$output" or die "COULD NOT OPEN $output :: $!\n";
}

if ($routput) {
	open ROUTPUT, ">$routput" or die "COULD NOT OPEN $routput :: $!\n";
}



$/ = ">";
eval ($input);
while (my $entry = eval ($input)){
	chomp($entry);
	#use only lower case letters to make it simple
	$entry =~ tr/ATCG/atcg/;
	(my $name) = split (/\n/, $entry);
	$entry =~ s/\n//g;
	$entry = substr ($entry, length($name));
	$entry = &clean ($entry);
	&process_concatamer ($minimum_size, $maximum_size, $entry, $mono);

}

my $entries = @sizes;
my $accumulated_total = 0;
&printerr ("MAX SIZE=$entries\n");
for (my $i = 0; $i < $entries; $i++){
	if ($sizes[$i] != 0){
		$accumulated_total += $sizes[$i];
		&printerr ("size $i: $sizes[$i] entries(up to now:$accumulated_total)\n");
	}
}


sub process_concatamer {
	my $minimum_size = shift;
	my $maximum_size = shift;
	my $concatamer = shift;
	my $mono = shift;
	&printerr ("===CONCATAMER:\n$concatamer\n");
	#to control the last entry
	$concatamer .= "X";
	my @all_ditags = split($punc,$concatamer);
	my $size = @all_ditags - 2;
	my $current = 0;
	my $ditag = shift (@all_ditags);
	my $separator = $punc;
	my $comment = $ditag;

	#
	#Processing the first ditag, may be a monotag
	#	
	if (length ($ditag) != 0) {
		if ($mono) {
			&process_monotag ($minimum_size, $maximum_size, $ditag, $mono, 0);
		}
		else {
			&printerr ("\tHEAD:$ditag\n");
		}
	}
	$ditag = shift (@all_ditags);
	while ($current < $size) {
		my $size_ditag = length($ditag);
		if (($minimum_size <= $size_ditag) && ($size_ditag <= $maximum_size)) {
			&printerr ("\tGOOD DITAG:$ditag\n");
			&printout ($ditag,"\n");
			$sizes[$size_ditag] += 1;
			$current += 1;
			$ditag = shift(@all_ditags);
			$comment = $ditag
		}
		elsif ($size_ditag > $maximum_size){
			if ($lost) {
				&process_long_ditags($minimum_size,$maximum_size,$ditag);
			}
			else {
				&printerr ("\tBAD DITAG(Large:$size_ditag):$ditag\n");
			}
			$current += 1;
			$ditag = shift(@all_ditags);
			$comment = $ditag
		}	 
		else {

			#
			#ditag is too small
			#
			my $next_ditag = shift(@all_ditags);
			my $size_next_ditag = length($next_ditag);
			my $combined_size = $size_ditag + 4 + $size_next_ditag;
			if ((($minimum_size/2) <= ($size_ditag+3)) && 
				(($size_ditag+1) <= ($maximum_size/2)) &&
				(($minimum_size/2) <= ($size_next_ditag+3)) && 
				(($size_next_ditag+1) <= ($maximum_size/2)) &&
				($minimum_size <= $combined_size) && 
				($combined_size <= $maximum_size)){
				#
				#seems that we have catg in the middle of a ditag
				#
				$ditag .= "$separator"."$next_ditag";
				$comment .= "[catg]$next_ditag";
				$sizes[$combined_size] += 1;
				&printerr ("\tJOINED DITAG:$ditag($comment)\n");
				&printout ("$ditag\n");
				#
				#get next ditag
				#
				$ditag = shift(@all_ditags);
				$current += 2;
				$comment = $ditag;

			}
			else {
				#
				#we can not accept catg in the middle otherwise
				# discard ditag
				&printerr ("\tBAD DITAG(Small:$size_ditag,$size_next_ditag):$ditag($comment)\n");
				$sizes[$size_ditag] += 1;
				$ditag = $next_ditag;
				$current += 1;
				$comment = $ditag;
			}
		}
	}
	my $last_size = length ($ditag);
	if ($last_size <= 1) {
		next;
	}
	$ditag = substr ($ditag, 0, $last_size - 1);
	
	if (length ($ditag) == 0) {
		return;
	}
	my $rest;
	if ($lost) {
		my $result = &process_long_ditags ($minimum_size,$maximum_size,$ditag);
		$rest = substr ($ditag, $result);
	}
	else {
		$rest = $ditag;
	}
	if ($mono) {
		&process_monotag ($minimum_size, $maximum_size, $rest, $mono, 1);
	}
	else {
		&printerr ("\tTAIL:$ditag\n");
	}
}



sub process_long_ditags {
	my $minimum_size = shift;
	my $maximum_size = shift;
	my $large_ditag = shift;
	my @spunc = split (//, $punc);
	my %ref = ("a", "[ctg]",
				"c", "[atg]",
				"t", "[acg]",
				"g", "[act]");
	
	my $variations = "(";
	my $punc_size = int (@spunc);
	$variations .= $ref{$spunc[0]} . substr ($punc, 1);
	for (my $i = 1; $i < $punc_size; $i++) {
		$variations .= "|" . substr ($punc, 0, $i) . $ref{$spunc[$i]} . substr ($punc, $i + 1);
	}
	$variations .= ")";
	my @all_fields = split(/$variations/,$large_ditag);
	my $size = @all_fields;
	my $current = 0;

	my $size = length ($large_ditag);
	&printerr ("\tLARGE DITAG($size):");
	my $sep = 0;
	foreach my $part (@all_fields) {
		if ($sep) {
			&printerr ("[$part]");
			$sep = 0;
		}
		else {
			&printerr ("$part");
			$sep = 1
		}
	}
	&printerr ("\n");

	my $ditag = shift(@all_fields);
	my $separator = shift(@all_fields);
	my $comment = $ditag;
	my $last_good = 0;

	LOOP:
	while ($current < $size){
		my $size_ditag = length($ditag);
		if (($minimum_size <= $size_ditag) && ($size_ditag <= $maximum_size)) {
			&printerr ("\t\tRECOVERED DITAG:$ditag($comment)\n");
			&printout ($ditag,"\n");
			$sizes[$size_ditag] += 1;
			$current += 2;
			$ditag = shift(@all_fields);
			$separator = shift(@all_fields);
			$comment = $ditag;
			$last_good += $size_ditag + 4;
		}
		elsif ($size_ditag > $maximum_size){
			#something is REALLY wrong!!
			&printerr ("\t\tBAD DITAG(Large:$size_ditag):$ditag($comment)\n");
			$sizes[$size_ditag] += 1;
			$current += 2;
			$ditag = shift(@all_fields);
			$separator = shift(@all_fields);
			$comment = $ditag;
		}	 
		else {

			#
			#ditag is too small
			#
			my $next_ditag = shift(@all_fields);
			my $next_separator = shift(@all_fields);
			my $size_next_ditag = length($next_ditag);
			#
			#separator is NOT "catg"
			#   in this case we keep joining ditags until correct size
			#
			$ditag .= $separator.$next_ditag;
			$comment .= "[$separator]$next_ditag";
			$separator = $next_separator;
			$current += 2;
		}
	}

	return $last_good;
}


sub process_monotag {
	my $minimum_size = shift;
	my $maximum_size = shift;
	my $entry = shift;
	my $mono = shift;
	my $side = shift;
	#0 for reverse, 1 to normal


	#
	#see if a ditag may be cutted of at the start
	#
	my $size_ditag = length ($entry);
	if ($size_ditag <= $maximum_size || $mono == 2) {
		#
		#check the size
		#
		my $monotag;
		if ($side) {
			$monotag = substr ($entry, 0, $maximum_size/2);
		}
		else {
			$monotag = substr ($entry, -$maximum_size/2);
		}
		if (length ($monotag) >= $minimum_size/2) {
			my @mask = ();
			for (my $i = 0; $i < $maximum_size/2; $i++) {
				push (@mask, "X");
			}
			#
			#put it on
			#
			&printerr ("\tMONOTAG:$monotag($entry)\n");
			if ($side) {
				&printout ($monotag, @mask, "\n");
			}
			else {
				&printout (@mask, $monotag, "\n");
			}
		}
		else {
			&printerr ("\tTOO SHORT MONOTAG:$entry\n");
		}
	}
	else {
		&printerr ("\tTOO LONG MONOTAG:$entry\n");
	}
}


sub printout {
	if ($output) {
		while (my $p = shift()) {
			print OUTPUT $p;
		}
	}
	else {
		while (my $p = shift()) {
			print $p;
		}
	}
}


sub printerr {
	if ($routput) {
		while (my $p = shift()) {
			print ROUTPUT $p;
		}
	}
	else {
		while (my $p = shift()) {
			print STDERR $p;
		}
	}
}


sub clean {
	my $string = shift;

	$string =~ s/^\s+//;
	$string =~ s/\s+$//;

	return $string;
}


