#!usr/bin/perl

use strict;
my %reads;
my $k = 0;

#Command line: <perl jxn_ID_1.0.pl (junction read bed file)>

open(BUFFER, $ARGV[0]);
while(my $line = <BUFFER>){
    my @fields = split("\t", $line);
#    $reads{$k}{name} = $fields[4];
    $reads{$k}{chr} = $fields[2];
    $reads{$k}{start} = $fields[3];
    $reads{$k}{CIGAR} = $fields[5];
    $k++;
}

close(BUFFER);

my %junctions;
my $j = 0;

LINE: foreach my $jxn(keys %reads){
    my @bases = split("M", $reads{$jxn}{CIGAR});
    my @bases2 = split("N", $bases[1]);
    my $five = $reads{$jxn}{start}+$bases[0];
    my $three = $reads{$jxn}{start}+$bases[0]+$bases2[0];
    foreach my $i(keys %junctions){
	if($reads{$jxn}{chr} == $junctions{$i}{chr}
	   && $five == $junctions{$i}{five}
	   && $three == $junctions{$i}{three}){
	    $junctions{$i}{count}++;
	    next LINE;
	}
    }
#    $junctions{$j}{name} = $reads{$jxn}{name};
    $junctions{$j}{chr} = $reads{$jxn}{chr};
    $junctions{$j}{five} = $five;
    $junctions{$j}{three} = $three;
    $junctions{$j}{count} = 1;
    $j++;
}

#my @counts;
#foreach my $m(keys %junctions){
#    unshift(@counts, $junctions{$m}{count}) unless $counts[0] > $junctions{$m}{count};
#}
#my $highcount = $counts[0];

open(OUT, ">>$ARGV[0].jxnquant.bed");
foreach my $l(sort{$junctions{$a}->{chr} cmp $junctions{$b}->{chr} ||
	$junctions{$a}->{five} <=> $junctions{$b}->{five} ||
	$junctions{$a}->{three} <=> $junctions{$b}->{three}} keys %junctions){
#    my $score = (1000*$junctions{$l}{count})/$highcount;
    print OUT "$junctions{$l}{chr}\t$junctions{$l}{five}\t$junctions{$l}{three}\t$junctions{$l}{count}\n";
}
close(OUT);
