#!/usr/bin/env perl
use warnings;
use strict;

use Getopt::Std;

# Convert variant calls in VCF format to MapSNPs input format.

my %opts;
getopts('gif:q:', \%opts) or die "ERROR: Wrong argument(s)\n";
my $genotype = $opts{'g'};
my $snpid    = $opts{'i'};
my $filter   = $opts{'f'} || ''; # filter string (exact match)
my $quality  = $opts{'q'} || 0;  # QC score threshold

my $line = 0; # input line counter
while (<>) {
  $line++;
  next if /^#/ || /^\s*$/;
  chomp;
  my @a = (split /\t/, $_, 11)[0..9]; # we only use 10 first columns
  @a = grep { defined } @a;
  my $colcount = scalar @a;
  unless ($colcount >= 8) {
    warn "WARNING: ($line) Insufficient number of fields ($colcount):\t" . join("\t", @a) . "\n";
    next;
  }
  map { s/^\s+//; s/\s+$//; s/^\.$// } @a;
  if ($genotype && !$a[9]) {
    warn "WARNING: ($line) Missing genotype field:\t" . join("\t", @a) . "\n";
    next;
  }
  # check for both alleles present
  warn("WARNING: ($line) Missing REF allele: " . join("\t", @a) . "\n"), next unless length $a[3];
  warn("WARNING: ($line) Missing ALT allele: " . join("\t", @a) . "\n"), next unless length $a[4];
  # skip indels and other possible multi-nucleotide variants
  next if length($a[3]) > 1;
  # validate REF and ALT alleles
  warn("WARNING: ($line) Invalid REF allele: " . join("\t", @a) . "\n"), next if $a[3] !~ /[ACTG]/i;
  warn("WARNING: ($line) Invalid ALT allele: " . join("\t", @a) . "\n"), next unless $a[4] =~ m|^[ACGT](?:[,][ACGT])*$|i;
  next if $filter && $filter ne $a[6];
  next if $quality && $a[5] < $quality;
  my $chrom = lc $a[0];
  $chrom = 'chr' . $chrom unless $chrom =~ /^chr/;
  $chrom =~ s/^chr(.+)$/chr\U$1/;
  print join("\t", $chrom . ':' . $a[1], $a[3] . '/' . $a[4]);
  my @comments;
  # genotype (first one only)
  push @comments, $a[9] if $genotype && length $a[9];
  # optional unique variant identifier, e.g. dbSNP rsID
  push @comments, $a[2] if $snpid && $a[2];
  print "\t# " . join("\t", @comments) if @comments;
  print "\n";
}
