#!/usr/local/bin/perl
# File: createSetBands
# Author: Terry Furey
# Date: 11/16/00
# Project: Human
# Description: Creates the set band file used for predicting chromomsome
#              band boundaries

# Usage message
if ($#ARGV != 3) {
  print stderr "USAGE: createSetBands <pct file> <insert file> <ctg dir> <window>\n";
  exit(1);
}

# Read parameters
$pctfile = shift(@ARGV);
$insfile = shift(@ARGV);
$ctgdir = shift(@ARGV);
$basewind = shift(@ARGV);
$window = $basewind * 1000;

# Determine the golden path positions for each of the inserts
open(INSERT, "<$insfile") || die("Could not open $insfile\n");
$line = <INSERT>; #header
while ($line = <INSERT>) {
  chomp($line);
  ($chr, $first, $second, $length, $type) = split(' ',$line);
  @allchrs = (@allchrs, $chr);
  @alltype = (@alltype, $type);
  @alllength = (@alllength, $length);  

  $chr = substr($chr,3);
  open(LIFT, "<$ctgdir/$chr/lift/ordered.lft") || die("Could not open $ctgdir/$chr/lift/ordered.lft\n");
  # Determine first window for the insert
  if ($first eq "-") {
    print stderr "Finding short arm chr$chr\n";
    @firstwind = (@firstwind, 0);
  } else {
    $found = 0;
    print stderr "Finding chr$chr $first\n";
    while(!$found) {
      $line = <LIFT>;
      chomp($line);
      ($ctgstart, $thisctg, $ctglen, $ctgchr, $chrlen) = split(' ',$line);
      ($dummy, $ctg) = split('/',$thisctg);
      if ($ctg eq $first) {
	$newwind = int(($ctgstart + $ctglen)/$window);
	@firstwind = (@firstwind, $newwind);
	$found = 1;
      }
    }
  }
  close(LIFT);
}
close(INSERT);

# Now, match up with PCT
open(PCT, "<$pctfile") || die("Could not open $pctfile\n");
while  ($line = <PCT>) {
  print stderr "Processing $line\n";
  chomp($line);
  ($chr, $type, $lines) = split(' ',$line);
  for ($i = 0; $i <= $#allchrs; $i++) {
    if (($allchrs[$i] eq $chr) && ($alltype[$i] eq $type)) {
      $index = $i;
      $i = $#allchrs + 1;
    }
  }
  if (($allchrs[$index] ne $chr) || ($alltype[$index] ne $type)) {
    die("Files not in sync, $chr $allchrs[$index] $type $alltype[$index]\n");
  }
  $numWindows = $alllength[$index]/$window;
  $bandstart = ($firstwind[$index] * $window);
  $lastend = $bandstart + $alllength[$index] - $window;
  for ($j = 0; $j < $lines; $j++) {
    $line = <PCT>;
    if (substr($line,0,1) ne "#") {
      chomp($line);
      ($chr, $band, $pctlen, $pct, $stain) = split(' ',$line);
      $bandwind = ($pct * ($alllength[$index]/$window));
      if (($bandwind%$basewind) < ($basewind/2)) {
	$bandwind -= $basewind;
      }
      if ($j == ($lines - 1)) {
	$bandend = $lastend;
      } else {
	$bandwind = int($bandwind/$basewind);
	$bandend = $bandstart + ($bandwind * $window);
      }
      print "$chr\t$band\t$bandstart\t$bandend\n";
      $bandstart = $bandend + $window;
    }
  }
}
close(PCT);
