#!/usr/bin/perl -w

die("Usage: mapview.pl <in.map>\n") if (@ARGV == 0);
my ($fh, $buf, @ref);
open($fh, "gzip -dc $ARGV[0] |") || die("[maq] Cannot decompress $ARGV[0]");
binmode($fh);
# read header
$buf = '';
read($fh, $buf, 8);
my ($format, $n_ref) = unpack("i2", $buf);
for (1..$n_ref) {
  read($fh, $buf, 4);
  my $len = unpack("i", $buf);
  read($fh, $buf, $len);
  ($_) = unpack("Z$len", $buf);
  push(@ref, $_);
}
read($fh, $buf, 8); # what is read is the number of reads, but frequently it is zero.
# read alignment
my @ACGT = ('A', 'C', 'G', 'T');
my $max = 127;
while (read($fh, $buf, 57+$max) == 57+$max) {
  my ($sq, $qse, $size, $q, $i1, $i2, $c0, $c1, $f, $aq, $sid, $x, $d, $rn)
	= unpack("a${max}cC8I2iZ36", $buf);
  print "$rn\t$ref[$sid]\t", ($x>>1)+1, "\t", ($x&1)?'-':'+', "\t$d\t$f\t$q\t$qse\t$aq\t",
	$i1&0xf, "\t$i2\t$c0\t$c1\t$size";
  my @t = unpack("C$size", $sq);
  my $seq = '';
  foreach my $p (@t) {
	if ($p == 0) { $seq .= 'N'; }
	else { $seq .= $ACGT[$p>>6]; }
	$p = chr(33 + ($p&0x3f));
  }
  print "\t$seq\t", @t, "\n";
}
close($fh);
