/* -*- mode: c -*-
 * $Id: moments.gcc,v 1.7 2001/01/05 02:07:51 doug Exp $
 * http://www.bagley.org/~doug/shootout/
 */

#include <stdio.h>
#include <stdlib.h>
#include <math.h>

#define MAXLINELEN 128

int
compare_doubles (const double *a, const double *b) {
    return (int) (*a - *b);
}


int
main() {
    char line[MAXLINELEN];
    int i, n = 0, mid = 0;
    double sum = 0.0;
    double mean = 0.0;
    double average_deviation = 0.0;
    double standard_deviation = 0.0;
    double variance = 0.0;
    double skew = 0.0;
    double kurtosis = 0.0;
    double median = 0.0;
    double deviation = 0.0;
    int array_size = 4096;

    double *nums = (double *)malloc(array_size * sizeof(double));

    while (fgets(line, MAXLINELEN, stdin)) {
    sum += (nums[n++] = atof(line));
    if (n == array_size) {
        array_size *= 2;
        nums = (double *)realloc(nums, array_size * sizeof(double));
    }
    }
    mean = sum/n;
    for (i=0; i<n; i++) {
    deviation = nums[i] - mean;
    average_deviation += fabs(deviation);
    variance += pow(deviation,2);
    skew += pow(deviation,3);
    kurtosis += pow(deviation,4);
    }
    average_deviation /= n;
    variance /= (n - 1);
    standard_deviation = sqrt(variance);
    if (variance) {
    skew /= (n * variance * standard_deviation);
    kurtosis = (kurtosis/(n * variance * variance)) - 3.0;
    }

    qsort(nums, n, sizeof (double),
      (int (*)(const void *, const void *))compare_doubles);
    mid = (n/2);
    median = n % 2 ? nums[mid] : (nums[mid] + nums[mid-1])/2;

    free(nums);

    printf("n:                  %d\n", n);
    printf("median:             %f\n", median);
    printf("mean:               %f\n", mean);
    printf("average_deviation:  %f\n", average_deviation);
    printf("standard_deviation: %f\n", standard_deviation);
    printf("variance:           %f\n", variance);
    printf("skew:               %f\n", skew);
    printf("kurtosis:           %f\n", kurtosis);
    return(0);
}