#!/usr/bin/perl

use Getopt::Std;
use strict;

my %opts = (p=>5, q=>0, c=>0, s=>undef);
getopts('p:c:sq:', \%opts);
my $p_num=$opts{p};
my $qu=$opts{q};
my $cov=$opts{c};
my $issdi=$opts{s};

die "Usage: plupsyn.pl [-p -c -s] Chromsfile pileupfile coorfile <chrid>\n" if (@ARGV <2); 
my $chromsfile = shift;
my $pileupfile = shift;
my $coorfile =shift;
print STDERR "reference chromosome file: ",$chromsfile, "\n";

my @chrs;
if (@ARGV) {
	@chrs=@ARGV;
} else {
	my @match_line;
###	fgrep_into { /^>/ } \@match_line, glob $chromsfile;
####	my @match_line=values %{$raw_fas[0]{matches}};
	open my $hfas, "grep '^>' $chromsfile |";
	while ( <$hfas> ) {
     		my $chid =~ s/^\s+//; 
     		($chid)=  $_ =~ m/^>(\S+)/;
     		push(@chrs, $chid);
	}
}	


foreach (@chrs) {
	my %seen;
	my $chid="^".$_.'[[:space:]]';
#	my $chid="^".$_;
	#my @raw_match;
	#fgrep_into { /$chid/ } \@raw_match, $pileupfile; 
	#foreach (@raw_match) {
	open my $hplup, "grep '$chid' $pileupfile |" ;
	while (<$hplup>) {
      		chomp;
      		my @F=split;
		if ($F[9] =~ /NOCONFIRM:([A-Z])/ ) {
		    $F[9]=$1;
		}
		my $sst;
		if ( (($F[3] =~/[Cc]/ && $F[9] eq "T") || ($F[3] =~/[Gg]/ && $F[9] eq "A")) && $F[8]>20 )  {
      			$seen{"$F[0]\t$F[1]\t$F[3]"}="\t$F[8]"."/"."$F[9]";
		} elsif ($F[3] eq $F[9]) {
			$seen{"$F[0]\t$F[1]\t$F[3]"}="\t$F[3]";
		} else {
			$seen{"$F[0]\t$F[1]\t$F[3]"}="\t?";
		}
	}

	my $hc;
	open $hc, "<$coorfile";
	while (<$hc>) {
		next if ( $_ !~ /$chid/ ); 
		chomp;
		my @F=split;
		if ($issdi) {
		    print "$F[0]\t$F[1]\t$F[2]\t$F[3]", $seen{"$F[0]\t$F[1]\t$F[3]"} ?$seen{"$F[0]\t$F[1]\t$F[3]"}:"\t?", "\n";
		} else {
		    print "$F[0]\t$F[1]\t$F[2]", $seen{"$F[0]\t$F[1]\t$F[2]"} ?$seen{"$F[0]\t$F[1]\t$F[2]"}:"\t?", "\n";
		}
	}
	close $hc;
}


