#!/usr/bin/perl
#
#    Usage: mean_stdev.pl files
#
#    Computes the mean and error of data which can be in several files.
#    Each data point is associated with a label (temperature say), so each
#    line has two entries,  y      T, and we want to average the y's
#    seperately for each T. The number of data points for each T can be
#    different.
#
#    The program uses associative arrays in which an element is specified by
#    any variable (in this case T), not necessarily an integer starting at 0.
#
while(<>)	# loop over all lines in all the files
{
    @line = split;     # split the line an put it into an array called line
    if($line[0] != "#")   # exclude comment lines
    {
	$T = $line[0]; # the first element is the temperature
	$count{$T}++;  # we count the number of data points at this temperature
	$y{$T}[$count{$T}] = $line[1]; # the second element is the data point, y
    }
}

print "     T    n_point     y_av       error \n";
foreach $T (sort {$a  <=> $b} (keys %count)) # loop over the temperatures
{
    $y_av = 0; $y_err = 0;        # Initialization
    $n = $count{$T};
    for ($i = 1; $i <= $n; $i++)  # loop over the data for a given temperature
    {
        $y_av += $y{$T}[$i];
        $y_err += $y{$T}[$i]**2;  # Note the square is "**" as in fortran
    }
    $y_av = &mean($y_av, $n);    # illustrates use of a subroutine in perl
    $y_err = &error($y_av, $y_err, $n);
    printf "%8.4f  %6d   %8.4f   %8.4f\n", $T, $n, $y_av, $y_err;
}

sub mean  #subroutine to compute the mean from the total and number of elements
{
    local($xav, $n) = @_;
    $xav / $n;
}

sub error #subroutine to compute the error from the mean and sum of squares
{
    local($xav, $xavsq, $n) = @_;
    $xavsq /= $n;
    sqrt(abs($xavsq - $xav**2) / ($n - 0.9999));
}
