#!/usr/bin/perl

open(IN, "/ddn/gs1/project/mousemeth/DNAseq/analysis/shared/reference_indel_files/indel_mm9.relative_to_BL6.bed"); @in = <IN>; close(IN);
open(OUT, ">indel.mm9_vs_BL6.bed");
for ($i=0; $i<$#in+1; $i++) {
  chomp $in[$i]; ($chr, $p0, $p2, $id) = split/\t/, $in[$i];
  if ($i == $#in) { $chrNext = ""; }
  else { chomp $in[$i+1]; ($chrNext, $p0Next, $p2Next, $idNext) = split/\t/, $in[$i+1]; }
  if (($p0 == $p0Next || $p2 == $p2Next) && $chr eq $chrNext) {
    if ($id eq "insertion" && $idNext =~ /deletion/) {
      @ar = split/\:/, $idNext; $ntI = commify($ar[1]);
      $delsize = $p2-$p0; $ntD = commify($delsize);
      print OUT "$chr\t$p0\t$p2\tBL6_indel_-$ntD/+$ntI\_nt\n"; 
    }
    elsif ($id =~ /deletion/ && $idNext eq "insertion") {
      @ar = split/\:/, $id; $ntI = commify($ar[1]);
      $delsize = $p2Next-$p0Next; $ntD = commify($delsize);
      print OUT "$chr\t$p0\t$p2Next\tBL6_indel_-$ntD/+$ntI\_nt\n";
    }
    else { print "ERROR: unexpected indel pair\n$in[$i]\n$in[$i+1]\n\n"; close(OUT); exit; }
    $i++;
  }

  else {
    if ($id eq "insertion") {
      $delsize = $p2-$p0; $nt = commify($delsize);
      print OUT "$chr\t$p0\t$p2\tBL6_deletion_-$nt\_nt\n";
    }
    else {
      @ar = split/\:/, $id; $nt = commify($ar[1]);
      $p0--;  # adjust for non-zero sized length
      print OUT "$chr\t$p0\t$p2\tBL6_insertion_+$nt\_nt\n";
    }
  }
}
close(OUT);
system "/ddn/gs1/home/grimmsa/tools/ucsc/bedToBigBed indel.mm9_vs_BL6.bed /ddn/gs1/project/mousemeth/reference/mm9/mm9assembly.chromSizes indel.mm9_vs_BL6.bigBed";


open(IN, "/ddn/gs1/project/mousemeth/DNAseq/analysis/shared/reference_indel_files/indel_mm9.relative_to_C3H.bed"); @in = <IN>; close(IN);
open(OUT, ">indel.mm9_vs_C3H.bed");
for ($i=0; $i<$#in+1; $i++) {
  chomp $in[$i]; ($chr, $p0, $p2, $id) = split/\t/, $in[$i];
  if ($i == $#in) { $chrNext = ""; }
  else { chomp $in[$i+1]; ($chrNext, $p0Next, $p2Next, $idNext) = split/\t/, $in[$i+1]; }
  if (($p0 == $p0Next || $p2 == $p2Next) && $chr eq $chrNext) {
    if ($id eq "insertion" && $idNext =~ /deletion/) {
      @ar = split/\:/, $idNext; $ntI = commify($ar[1]);
      $delsize = $p2-$p0; $ntD = commify($delsize);
      print OUT "$chr\t$p0\t$p2\tC3H_indel_-$ntD/+$ntI\_nt\n";
    }
    elsif ($id =~ /deletion/ && $idNext eq "insertion") {
      @ar = split/\:/, $id; $ntI = commify($ar[1]);
      $delsize = $p2Next-$p0Next; $ntD = commify($delsize);
      print OUT "$chr\t$p0\t$p2Next\tC3H_indel_-$ntD/+$ntI\_nt\n";
    }
    else { print "ERROR: unexpected indel pair\n$in[$i]\n$in[$i+1]\n\n"; close(OUT); exit; }
    $i++;
  }
  else {
    if ($id eq "insertion") {
      $delsize = $p2-$p0; $nt = commify($delsize);
      print OUT "$chr\t$p0\t$p2\tC3H_deletion_-$nt\_nt\n";
    }
    else {
      @ar = split/\:/, $id; $nt = commify($ar[1]);
      $p0--;  # adjust for non-zero sized length
      print OUT "$chr\t$p0\t$p2\tC3H_insertion_+$nt\_nt\n";
    }
  }
}
close(OUT);
system "/ddn/gs1/home/grimmsa/tools/ucsc/bedToBigBed indel.mm9_vs_C3H.bed /ddn/gs1/project/mousemeth/reference/mm9/mm9assembly.chromSizes indel.mm9_vs_C3H.bigBed";









sub commify {
  $number = $_[0];
  $number =~ s/(\d)(?=(\d{3})+(\D|$))/$1\,/g;
  return $number;
}