#!/usr/bin/perl

use strict;


my @chrs = (1..19);
print("chrs = @chrs\n");

my ($line, $start, $ch);

my $ofile = 'count.gatk.raw.txt';

open(OUT, '>', $ofile) or die "Cannot open output file, reason $!\n";
print(OUT "Chr Total Type 101 211 218 227 231 234 236 246 303 310 311 313 314 327 329 330 334 34 44 54\n");

foreach $ch ( @chrs ) {
    print("chromosome $ch\n");

    my $pname = 'gatk/raw/snps.raw.chr'.$ch.'.vcf';
#    my $pname = 'gatk/raw/filter/output/snps.filtered.chr'.$ch.'.vcf';

    open(IN, '<', $pname) or die "Cannot open input file $pname, reason $!\n";

    $line = <IN>;
    $start = substr( $line, 0, 2);

    while ( $start eq '##' ) {
        $line = <IN>;
        $start = substr( $line, 0, 2);
    } # end of headers while loop
 
    my $header_line = $line;


    my ($i, $mouse, $pval);
    my (@row, @mouse_data);

    my $Ntot = 0;

    my @Nmiss = (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0);
    my @Nref = (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0);
    my @Nhet = (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0);
    my @Nalt = (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0);


    while ( defined($line = <IN>) ) {

	chomp($line);
	@row = split("\t", $line);

	if ( length($row[4]) == 1 ) {

	  $Ntot++;  

	  for ($i=0; $i<20; $i++) {

	    $mouse = $row[9+$i];
	    @mouse_data = split(':', $mouse);

	    my $geno = $mouse_data[0];

	    if ( $geno eq './.' ) {
	      $Nmiss[$i]++;
	    } elsif ( $geno eq '0/0' ) {
	      $Nref[$i]++;
	    } elsif ( $geno eq '1/1' ) {
	      $Nalt[$i]++;
	    } elsif ( $geno eq '1/0' || $geno eq '0/1' ) {
	      $Nhet[$i]++;
	    } else {
	      print("ERROR - unrecognised genotype geno $geno  bp $row[1]  mouse $i  data $mouse\nrow @row\nline $line\n");
	      exit;
	    }

	  } # end of loop over mice

	} # end of mult ALT loop

    } # end of while loop
    close(IN);

    print(OUT "$ch $Ntot MISS @Nmiss\n");
    print(OUT "$ch $Ntot REFS @Nref\n");
    print(OUT "$ch $Ntot HETS @Nhet\n");
    print(OUT "$ch $Ntot ALTS @Nalt\n");

} # end of chr loop

close(OUT);
