/* mean_stddev 
   -- computes the mean and variance of an array a[n]
   -- returns mean, sets variance argument
   -- uses standard one-pass algorithm */
float mean_var(float a[], int n, float *variance)
{
  int  i;
  float mean, sum, sum_squares;

  sum = sum_squares = 0.0;
  for ( i = 0; i < n; i++ )
    {
      sum += a[i];
      sum_squares += a[i]*a[i];
    }

  mean = *variance = 0.0;
  if ( n > 0 )
    mean = sum / n;
  if ( n > 1 )
    *variance = (sum_squares - n*mean*mean) / (n-1);

  return mean;
}

/* mean_var2
   -- computes the mean and variance of an array a[n]
   -- returns mean, sets variance argument
   -- uses standard two-pass algorithm */
float mean_var2(float a[], int n, float *variance)
{
  int  i;
  float mean, sum, sum_squares;

  sum = sum_squares = 0.0;
  for ( i = 0; i < n; i++ )
    sum += a[i];

  if ( n > 0 )
    mean = sum / n;

  if ( n > 1 )
    {
      for ( i = 0; i < n; i++ )
        sum_squares += (a[i]-mean)*(a[i]-mean);
      *variance = sum_squares / (n-1);
    }

  return mean;
}