|
Packit |
d3f73b |
/*
|
|
Packit |
d3f73b |
* Experimental data distribution table generator
|
|
Packit |
d3f73b |
* Taken from the uncopyrighted NISTnet code (public domain).
|
|
Packit |
d3f73b |
*
|
|
Packit |
d3f73b |
* Rread in a series of "random" data values, either
|
|
Packit |
d3f73b |
* experimentally or generated from some probability distribution.
|
|
Packit |
d3f73b |
* From this, report statistics.
|
|
Packit |
d3f73b |
*/
|
|
Packit |
d3f73b |
|
|
Packit |
d3f73b |
#include <stdio.h>
|
|
Packit |
d3f73b |
#include <stdlib.h>
|
|
Packit |
d3f73b |
#include <math.h>
|
|
Packit |
d3f73b |
#include <malloc.h>
|
|
Packit |
d3f73b |
#include <sys/types.h>
|
|
Packit |
d3f73b |
#include <sys/stat.h>
|
|
Packit |
d3f73b |
|
|
Packit |
d3f73b |
void
|
|
Packit |
d3f73b |
stats(FILE *fp)
|
|
Packit |
d3f73b |
{
|
|
Packit |
d3f73b |
struct stat info;
|
|
Packit |
d3f73b |
double *x;
|
|
Packit |
d3f73b |
int limit;
|
|
Packit |
d3f73b |
int n=0, i;
|
|
Packit |
d3f73b |
double mu=0.0, sigma=0.0, sumsquare=0.0, sum=0.0, top=0.0, rho=0.0;
|
|
Packit |
d3f73b |
double sigma2=0.0;
|
|
Packit |
d3f73b |
|
|
Packit |
d3f73b |
fstat(fileno(fp), &info;;
|
|
Packit |
d3f73b |
if (info.st_size > 0) {
|
|
Packit |
d3f73b |
limit = 2*info.st_size/sizeof(double); /* @@ approximate */
|
|
Packit |
d3f73b |
} else {
|
|
Packit |
d3f73b |
limit = 10000;
|
|
Packit |
d3f73b |
}
|
|
Packit |
d3f73b |
x = (double *)malloc(limit*sizeof(double));
|
|
Packit |
d3f73b |
|
|
Packit |
d3f73b |
for (i=0; i
|
|
Packit |
d3f73b |
fscanf(fp, "%lf", &x[i]);
|
|
Packit |
d3f73b |
if (feof(fp))
|
|
Packit |
d3f73b |
break;
|
|
Packit |
d3f73b |
sumsquare += x[i]*x[i];
|
|
Packit |
d3f73b |
sum += x[i];
|
|
Packit |
d3f73b |
++n;
|
|
Packit |
d3f73b |
}
|
|
Packit |
d3f73b |
mu = sum/(double)n;
|
|
Packit |
d3f73b |
sigma = sqrt((sumsquare - (double)n*mu*mu)/(double)(n-1));
|
|
Packit |
d3f73b |
|
|
Packit |
d3f73b |
for (i=1; i < n; ++i){
|
|
Packit |
d3f73b |
top += ((double)x[i]-mu)*((double)x[i-1]-mu);
|
|
Packit |
d3f73b |
sigma2 += ((double)x[i-1] - mu)*((double)x[i-1] - mu);
|
|
Packit |
d3f73b |
|
|
Packit |
d3f73b |
}
|
|
Packit |
d3f73b |
rho = top/sigma2;
|
|
Packit |
d3f73b |
|
|
Packit |
d3f73b |
printf("mu = %12.6f\n", mu);
|
|
Packit |
d3f73b |
printf("sigma = %12.6f\n", sigma);
|
|
Packit |
d3f73b |
printf("rho = %12.6f\n", rho);
|
|
Packit |
d3f73b |
/*printf("sigma2 = %10.4f\n", sqrt(sigma2/(double)(n-1)));*/
|
|
Packit |
d3f73b |
/*printf("correlation rho = %10.6f\n", top/((double)(n-1)*sigma*sigma));*/
|
|
Packit |
d3f73b |
}
|
|
Packit |
d3f73b |
|
|
Packit |
d3f73b |
|
|
Packit |
d3f73b |
int
|
|
Packit |
d3f73b |
main(int argc, char **argv)
|
|
Packit |
d3f73b |
{
|
|
Packit |
d3f73b |
FILE *fp;
|
|
Packit |
d3f73b |
|
|
Packit |
d3f73b |
if (argc > 1) {
|
|
Packit |
d3f73b |
fp = fopen(argv[1], "r");
|
|
Packit |
d3f73b |
if (!fp) {
|
|
Packit |
d3f73b |
perror(argv[1]);
|
|
Packit |
d3f73b |
exit(1);
|
|
Packit |
d3f73b |
}
|
|
Packit |
d3f73b |
} else {
|
|
Packit |
d3f73b |
fp = stdin;
|
|
Packit |
d3f73b |
}
|
|
Packit |
d3f73b |
stats(fp);
|
|
Packit |
d3f73b |
return 0;
|
|
Packit |
d3f73b |
}
|