#!/usr/bin/perl
#
#This script accepts 3 arguments, minimum_size, maximum_size and file
#the first two arguments stat the minimum and maximum ditag size (considering that
#tags can have different sizes. 
#

use strict;

my $min_tag_size = 17;
my $max_tag_size = 18;
my $coesive = 0;
my $input = "";
my $output = "out";
my $routput = "report";

my $version = 0.43;
my $vname = "sage_analysis.pl";

if (!$ARGV[0] || $ARGV[0]eq "-h") {
	print "\n\t-min <integer> \t:: set minimum size of the tags [default:$min_tag_size]\n",
	"\t-max <integer> \t:: set maximum size of the tags [default:$max_tag_size]\n",
	"\t-c <integer> \t:: set # of overlapping bases between tags [default:$coesive]\n",
	"\t-i <FILE> \t:: input file name [default:STDIN]\n",
	"\t-o <FILE> \t:: output file name prefix [default:$output]\n",
	"\t-r <FILE> \t:: report file name prefix [default:$routput]\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)) {
	#
	#coesive
	#
	if ($args eq "-c") {
		if ($args = shift(@ARGV)) {
			$coesive = $args;
		}
		else {
			$coesive = 2;
		}
	}
	#
	#minimun size
	#
	elsif ($args eq "-min") {
		if ($args = shift(@ARGV)) {
			$min_tag_size = $args;
		}
		else {
			print STDERR "MINIMUM SIZE NOT SPECIFIED!\n";
			exit;
		}
	}

	#	
	#maximun size
	#
	elsif ($args eq "-max") {
		if ($args = shift(@ARGV)) {
			$max_tag_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; 
		}
	}


}

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

#
#compute limits for valid ditag size
#
my $start = 2 * $min_tag_size;
my $end = 2 * $max_tag_size ;

if ($coesive){
	$start -= $coesive;
	$end -= $coesive;
}

#
#first generate ditags and statitistics on concatamers
#
print STDERR ("extract_ditags.pl -min $start -max $end -i $inputO -o $output\_ditags -r $routput\_ditags\n");
system("extract_ditags.pl $start $end -i $inputO -o $output\_ditags -r $routput\_ditags");

#
#now compute ditags
#
print STDERR ("extract_tags.pl -min $min_tag_size -max $max_tag_size -c $coesive -i $output\_ditags -o $output\_tags -r $routput\_tagsStat$start$end\n");
system("extract_tags.pl -min $min_tag_size -max $max_tag_size -c $coesive -i $output\_ditags -o $output\_tags -r $routput\_tags");

#
#if needed, print statistics
#

my $complete_report = 0;
if ($complete_report) {
	my $ln;
	my $trash;
	system("grep GOOD erro$start$end | wc > temp");
	system("grep JOIN erro$start$end | wc >> temp");
	system("grep BAD erro$start$end | grep Small |wc >> temp");
	system("grep BAD erro$start$end | grep Large |wc >> temp");
	open(TEMP,"temp");
	$ln = <TEMP>;
	($trash, my $good) = split(/[\s\t]+/,$ln);
	$ln = <TEMP>;
	($trash, my $join) = split(/[\s\t]+/,$ln);
	$ln = <TEMP>;
	($trash, my $small) = split(/[\s\t]+/,$ln);
	$ln = <TEMP>;
	($trash, my $large) = split(/[\s\t]+/,$ln);
	close(TEMP);
	open (TEMP,"tagsStat$start$end");
	$/= "====\n";
	my $data = <TEMP>;
	chomp($data);
	close(TEMP);
	$/= "\n";
	my @dataLines = split("\n",$data);
	my $lines = @dataLines;
	($trash,my $unique) = split("=",$dataLines[$lines-2]);
	($trash,my $conflicts) = split("=",$dataLines[$lines-1]);
	(my $keysInConflict,my $tagsInConflict,my $totalTags, my $totalGoodTags) = split(",",$conflicts);
	&printout ("$start"."x"."$end\t$good\t$join\t$small\t$large\t$unique\t$keysInConflict\t\t$tagsInConflict\t\t$totalTags\t\t$totalGoodTags\n");
	&printerr ("good=$good;bad=0;join=$join\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;
		}
	}
}
