#!/usr/bin/perl

# Time-stamp: <2005-07-07 14:27:16 kasper>

use warnings;
use strict;
use Getopt::Long;
use constant USAGE =><<END;

SYNOPSIS:

 tablestat.pl [OPTIONS] [infile [outfile]]

DESCRIPTION:

This script calculates average etc for numbers in column <col> lines
starting with '#' are ignored


OPTIONS:

       --help
             Prints this help.
       --column <int>
             Column to print stats for. Default is 1.

EXAMPLES:

 tablestat.pl --column 4 infile output

 cat infile | tablestat.pl > output

DEPENDENCIES

 None.

AUTHOR:

Anders Krogh / Kasper Munch

COPYRIGHT:

This program is free software. You may copy, modify, and redistribute
it under the same terms as Perl itself.

END

my $help = 0;
my $column = 1;

GetOptions(
           "help" => \$help,
           "column=i" => \$column,
          ) or die USAGE;

$help and die USAGE;

@ARGV = ('-') unless @ARGV;
my $input = shift @ARGV;
open my $in, "$input" or die "$input: $!\n";
@ARGV = ('>&STDOUT') unless @ARGV;
my $output = shift @ARGV;
open my $out, ">$output" or die "$output: $!\n";

my $C = $column - 1;

# Print an error message and exit
sub Error {
    print "$_[0]\n";
    exit(1);
}

my $min =  10000000000.;
my $max = -10000000000.;

my $sum = 0.;
my $ssum = 0.;
my $N = 0;

while (<$in>) {
    chop;
    if (/^[^#]/) {
	$_ =~ s/^\s+//;
        my @Fld = split(/\s+/, $_, 9999);
	next if ( $#Fld< $C );
	my $x = ( $Fld[$C] );
	$sum += $x;
	$ssum += $x*$x;
	$N += 1;
	$min = $x if ($x < $min);
	$max = $x if ($x > $max);
    }
}

my $av = $sum/$N;

my $var = $ssum/($N-1) - $N * $av * $av/($N-1);

my $stddev = sqrt($var);

print $out "Points:   $N\n";
print $out "Mean:     $av\n";
print $out "Variance: $var\n";
print $out "Stddev:   $stddev\n";
print $out "Max:      $max\n";
print $out "Min:      $min\n";
print $out "Sum:      $sum\n";


=head1 SYNOPSIS:

 tablestat.pl [OPTIONS] [infile [outfile]]

=head1 DESCRIPTION:

This script calculates average etc for numbers in column <col> lines
starting with '#' are ignored


=head1 OPTIONS:

=over 4

=item --help

Prints this help.

=item --column <int>

Column to print stats for. Default is 1.

=back

=head1 EXAMPLES:

 tablestat.pl --column 4 infile output

 cat infile | tablestat.pl > output

DEPENDENCIES

 None.

=head1 AUTHOR:

Anders Krogh / Kasper Munch

=head1 COPYRIGHT:

This program is free software. You may copy, modify, and redistribute
it under the same terms as Perl itself.


=cut